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PREFACE 


This  report  presents  the  work  since  the  previous  semiannual 
report  cn  Contract  Nonr  4647(00).  It  also  summarizes  the  previous  work 
on  the  contract  toward  the  development  of  quantitative  models  for  high 
energy  pulsed  arcs  of  the  type  that  are  used  for  the  optical  pumping 
of  high  energy  laiers.  The  studies  in  this  contract  have  been  concerned 
primarily  with  xenon  as  this  gas  has  been  used  most  successfully  for 
high  energy  laser  pumping;  the  principles  contained  in  the  models  for 
the  arcs  are  applicable  to  other  highly  radiative  arcs,  though  in  many 
cases  the  necessary  quantitative  descriptions  of  the  physical  properties 
of  the  arc  plasma  are  not  available  at  the  present  time. 

Recent  developments  in  pumping  geometries  and  in  additive 
lamps,  arising,  in  part,  from  this  work  and  from  work  in  a  closely 
related  area,  the  optical  pumping  of  continuous  lasers,  under  Contracts 
DA-28-043-AMC-00292(E)  and  DA-28-04 3-AMC-0209 7 (E) ,  indicate  new  promise 
for  efficient  optically  pumped  lasers.  The  principles  discussed  in  this 
report  provide  a  basis  for  an  understanding  of  the  optical  pumps  for 
lasers . 


ABSTRACT 


This  report  summarizes  the  studies  made  on  Contract  Nonr 
4647(00)  towards  the  development  of  models  for  the  highly  radiative 
arcs  used  for  the  high  energy  pumping  of  lasers.  The  report  alsi 
presents  the  experimental  and  theoretical  studies  since  the  last  semi¬ 
annual  report.  The  experimental  investigations  were  primarily  concerned 
with  more  extensive  measurements  of  the  speccral  radiance  of  the  plasma 
to  provide  verification  for  the  models.  The  theoretical  work 
has  resulted  in  computer  methods,  described  in  the  appendices,  to  calcu¬ 
late  the  transport  properties,  the  spectral  absorptivities  for  the 
lines  and  the  continuum  of  xenon,  and  the  spectral  radiance  and 
temperature  profiles  in  cylindrical  arcs.  Also  included  as  an  appendix 
is  a  theoretical  analysis  of  the  xenon  arc  using  radiative  transport 
techniques  developed  in  other  studies. 
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I.  INTRODUCTION 

The  optically  pumped  solid  state  laser  provides  the  highest 
peak  powe*-  and  highest  peak  energy  available  from  coherent  sources. 

This  report  summarizes  the  work  on  Contract  Nonr  4647(00)  toward  the 
development  of  quantitative  models  for  the  highly  radiative  arcs  used  to 
excite  these  lasers.  Contained  in  the  appendices  to  this  report  are 
computer  programs  developed  in  this  work  for  the  calculation  of  many 
features  of  these  arcs.  The  calculation  procedures  are  for  xenon, 
thou  h  in  principle,  they  are  applicable  to  other  atomic  systems. 

The  theoretical  and  the  experimental  research  on  this  contract 
has  been  directed  towards  the  creation  of  quantitative  models  for  the 
flash  iamps  that  in  turn  wo”ld  be  useful  in  the  improvement  of  high 
energy  lasers.  The  development  of  the  models  required  basic  studies  of 
the  many  aspects  of  highly  radiative  arcs,  and  lead  to  considerabl  j 
expansion  of  the  experimental  techniques  available  for  the  study  cf  these 
dense  plasmas.  These  plasmas  may  be  encountered  in  many  situations  other 
than  the  flash  tubes  for  pumping  high  energy  lasers,  including  re-entry 
shock  waves  and  high  voltage  switch  gear. 

The  effort  in  this  contract  can  be  divided  into  three  areas: 

(1)  the  development  of  theoretical  techniques  to  calculate  the  necessary 
physical  properties  of  arc  plasmas,  (2)  the  applications  of  these  tech¬ 
niques  to  create  models  for  the  high  radiative  arcs,  and  (3)  the  experi¬ 
mental  measurements  required  to  formulate  and  then  to  verify  the  models. 
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The  earlier  reports  (references  1,  2,  3  &  4)  described  many  aspects  of 
the  problem  more  fully.  The  computer  programs  resulting  from  this  work 
are  all  in  the  ALGOL  language  for  a  Burroughs  B5500  DISK  computer. 

The  development  of  the  models  required  techniques  for  cal¬ 
culating  the  transport  properties  (electrical  and  thermal  conductivity) , 
and  the  spectral  absorptivity  (the  broadened  spectral  lines  plus  the 
continuum)  as  a  function  of  temperature  and  pressure.  In  addition, 
methods  for  handling  energy  transport  within  highly  radiative  inhomo¬ 
geneous  temperature  nrc  plasmas  were  required. 

For  the  transport  properties  of  xenon,  we  were  fortunate  in 
obtaining  the  assistance  of  Dr.  R.  S.  de  Voto  of  Stanford  University, 
who  has  developed  computer  programs  for  calculating  the  properties  of 
partially  ionized  monatomic  gases.  His  program  for  the  calculation  of 
the  electrical  and  thermal  conductivities  of  xenon  as  a  function  of 
temperature  and  pressure  is  in  Appendix  A. 

The  calculation  of  th«.  spectral  absorptivity  of  xenon  required 
consideration  of  the  broadened  lines  and  the  continuum.  The  final  pro¬ 
gram  given  in  Appendix  B  is  due  to  the  efforts  of  many  people,  partic¬ 
ularly  D.  Sun  of  California  Institute  of  Technology,  E.  Corinaldesi,  now 
at  Boston  University,  and  R.  G.  Gchlecht,  now  at  Aeroneutronics  Division 
of  Ford  Motor  Company.  Appendix  C  uses  the  bound-bound  transition  prob¬ 
abilities  and  other  information  on  xenon  to  calculate  the  optical 
radiation  emitted  by  a  cylindrical  arc  of  a  given  temperature  profile 
and  pressure.  Appendix  D  presents  calculations  for  Whittaker  functions 
of  real  and  imaginary  argument  for  use  in  calculations  of  the  matrix 
elements  for  bound-free  and  free-free  transitions.  We  will  discuss 
details  of  these  calculations  in  a  later  section. 

Two  approaches  were  followed  in  the  creation  of  the  models  for 
the  arc.  In  one  approach,  we  assumed  the  arc  to  be  homogeneous  in  tem¬ 
perature  and  neglected  thermal  conduction  to  the  walls.  Experimental 
measurements  had  indicated  that  this  was  a  reasonable  assumption  for  the 
central  core.  Such  model  would  explain  a  number  of  the  observed  cha^ac- 
eristics,  such  as  the  arc  radiant  emission  saturation  in  the  infrared 
and  in  the  visible,  and  account  for  a  major  portion  of  the  input  power. 
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The  pressure  to  be  used  for  this  homogeneous  temperature  model  is  open 
to  question  due  to  the  actual  temperature  gradient  which  exists  near 
the  walls.  j.n  the  other  approach,  the  electric  field  across  the  arc,  the 
inside  diameter  of  the  arc  .ube,  and  the  pressure  of  heavy  particle 
number  density  are  to  be  specified  and  the  radial  absolute  temperature 
distribution  is  to  be  sought.  With  the  temperature  distribution,  we 
could  calculate  all  of  the  observed  properties  of  the  arc  including 
the  spectral  radiance  at  the  walls  (Appendix  C) .  This  calculation  re¬ 
required  the  solution  to  the  energy  balance  equation  for  a  non-grav, 
nonhomogeneous  temperature  plasma.  An  approximate  analysis  described 
in  Appendix  E  was  applied  to  this  problem  to  obtain  a  better  initial 
temperature  distribution.  The  final  program  developed  in  this  work  is 
given  in  Appendix  F. 

The  experimental  measurements  were  intended  to  aid  in  the 
formulation  and  verification  of  the  models.  In  the  course  of  the  work, 
the  flash  lamp  or  more  completely,  the  confined  wall-stabilized  pulsed 
arc  discharge,  was  found  to  be  a  useful  vehicle  to  study  high  power 
density,  high  pressure,  highly  radiative  plasmas.  Plasmas  with  these 
characteristics  In  gases  other  than  xenon  are  of  interest  to  workers 
in  many  fields  or  research  and  applications.  In  the  course  of  this  work 
and  in  the  closely  associated  company  sponsored  research,  techniques  have 
been  and  are  being  developed  to  study  these  plasmas.  These  techniques 
have  included  two  rapid  scanning  spectrometers,  one  of  high  resolution 
scanning  over  a  limited  range^  the  other  of  decreased  resolution  but 
with  a  much  wider  spectral  range^. 

A  problem  still  exis  *  in  the  accurate  measurement  of  the 
instantaneous  pressure.  The  pressure  is  crucial  to  an  accurate  descrip¬ 
tion  of  the  arc  as  the  pressure  varies  widely  with  the  temperature 
gradient  near  the  wall.  Currently  a  laser  interferometer  similar  to  that 

(o\ 

of  Gerardo  et  al.  is  being  set  up  to  allow  electron  density  measure¬ 
ments.  These  electron  density  measurements  will  serve  to  check  the 
values  we  have  assumed  in  the  calculations. 
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II.  BASI.  PRINCIPLES  OF  HIGHLY  RADIATIVE  ARCS 

The  arcs  used  for  pumping  high  energy  lasers  have  many 
features  in  common  with  arcs  used  for  continuous  laser  pumping.  The 
. .  B'  .ry  difference  between  these  arcs  lies  in  the  size  of  tne  arc 
column  and  the  power  density  in  the  arc.  In  both  cases,  the  arcs  are 
near  or  above  atmospheric  pressure.  As  the  efficiency  of  the  laser 
systems  increase,  the  differences  between  the  two  arc  applications 
are  becoming  less  marked. 

The  usual  corm  of  the  arc  used  in  high  energy  laser  pumping 
is  a  cylindrical  arc  surrounded  by  a  wall.  We  have  only  investigated 
this  wall  stabilized  arc.  In  the  arc,  current  flow  carries  power  into 
the  arc;  radiation  external  to  the  arc  and  thermal  conduction  to  the 
arc  surroundings  (i.e.:  the  wall)  carry  power  a vrv  from  the  arc. 
Convective  power  transport  can  be  neglected  in  this  type  of  wall  stabi¬ 
lized  confined  arc  unless  there  is  gas  flo” — which  we  have  not  considered 
in  this  work.  A  wal  stabilized  highly  radiative  arc  such  as  those  used 
for  high  energy  laser  pumping  has  a  central  plasma  core  at  the  tempera¬ 
ture  of  the  fully  ionized  plasma  (the  actual  temperature  varying  with 
the  composition  of  the  plasma)  with  the  temperature  at  the  wall  dependent 
upon  che  rate  at  which  power  is  carried  away  from  the  walls.  (A  fully 
ionized  plasma  may  be  only  .01%  ionized.  It  is  characterized  by  the 
doninance  of  the  electron-electron  and  electron-ion  interaction  over 
the  other  particle  interaction). 

The  relation:  between  the  input  pewer  and  that  carried  away 
from  the  arc  by  radi  ition  and  conduction  can  be  put  into  equation  form. 
Figure  1  presents  the  family  of  equations  that  relate  the  power  input 
per  unit  volume  to  the  powe’-  carried  away  by  radiation  and  -onduction. 

«e  consider  the  arc  to  be  in  local  thermal  equilibrium  (LTE)  and 
thus  descr^bable  by  Boltzmann's  equation,  Saha's  equation  and  a 
simple  radiation  law.  The  electrical  and  thermal  conductivities  as 
well  as  the  spectral  absorptivities  are  functions  of  the  local  tempera¬ 
ture  and  pre^  ure  for  each  constituent.  The  various  constituents 
of  the  plasma  interact  with  e*ch  other  thro:-  h  the  ionization  law 
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The  Energy  Balance 


(frad  +  ftc)  ■  c  (p'T)  E 


where  E  is  the  electric  field  in  Volts  cm  * ,  o  is  the  electrical  conductivity  in 


c-1 


(ohm  cm)  and  is  depsndent  uoon  the  pressure  (p)  and  the  temperature  (TI .  F. 

_2 

is  the  thermally  conducted  power  flux  in  Watts  cm  and  is  given  bv 


T.C. 


T.C. 


-K(p,T)  grad  T 


-2  -1 

K  is  the  thermal  conductivity  in  Watts  cm  °K  cm. 

_a 

(l.e.  radiative  flux)  in  Watts  cm  '  and  is  expressed  in  terms  of  F 
The  spectral  radiant  emittance 
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where  v  is  the  frequency  and  )  the  wavelength  in  consistent  mits.  Consider  the 

•  component  of  F  ;  it  is  expressed  in  terms  of  I  - 
VRAD  V 

F  -  j  1  (S)  Cos  (x,S)  d. 

'x  u 

I  is  the  spectral  radiance  (i.e.  specific  intensity)  in  Watts  cm~Z  Steradian 
V  -1 

frequency  unit  along  the  S  direction.  u  is  the  solid  angle  over  which  the 
Integral  is  to  be  taken.  1  is  the  solution  to  the  equation  of  transfer  for  the 
geometry  being  considered: 


dys) 

ds 


e  ' (d.T)  (B  (T)  -  I  (S)  ] 


in  which  s  is  the  magnitude  of  the  direction  vector  S  a. id  <  '  is  the  spectral 

-1 

absorptivity  including  stimulated  emission  in  cm  and  is  given  by 
i  i  -wi  -hv/»-T 

<  ■  x  (p, i) (1-e  ) 

V  V 

B  (T)  is  the  Plane.:  distribution  function,  h.  is  Planck's  constant,  *  is  Boltzmann's 

constant  and  are  to  be  in  units  consistent  with  T  and  v.  <  is  the  spectral 

-1  v 
absorptivity  in  cm  which  in  turn  is  the  sum  of  the  spectral  absorptiviti es 

for  bound-bound  transitions,  ^  ,  the  spectral  absorptivity  for  bound-free 

transitions,  <  and  the  spectral  absorptivity  for  free-free  transitions,* 

v  v 


* 

V 
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The  relations  between  the  spectral  absorptivi tf es  and  the  transition  probabilities 
for  these  processes  will  be  presentee,  in  Appendix  C. 


FIGURE  1 

The  equations  relating  the  electrical  power 
input  to  the  optical  radiation  and  thermal  conduction. 


(Saha  s  equation)  and  thus  the  electron  density,  and  in  broadening 
the  lines  through  electron,  ion,  and  Van  der  Waals  mechanisms. 

2 . 1  Electrical  and  Thermal  Conductivity 

The  transport  properties  (i.e.,  the  electrical  and  thermal 
conductivities)  of  partially  ionized  gases  are  currently  under  extensive 
investigation  in  many  laboratories.  We  have  used  de  Voto's  calcula¬ 
tions^’  as  applied  to  xenon  using  the  momentum  transfer  cro°s  sec¬ 
tion  of  Frost  and  Phelps.  The  program  for  the  electrical  and  ther¬ 

mal  conductivities  of  xenon  is  given  in  Appendix  A. 

2  2  The  Spectral  Absorptivity 

The  spectral  absorptivity  of  atomic  species  can  now  be  calcu¬ 
lated  (in  theory)  to  a  relatively  high  precision  without  exorbitant 
computer  time  using  "quantum  defect"  methods  to  obtain  the  transition 
probabilities  for  bound-bound  transitions  and  the  cross  sections  of 

the  bound-free  and  free-free  transitions.  The  electron  broadening  was 

(12) 

calculated  using  the  line  broadening  theory  outlined  by  Griem  as 

(13) 

applied  by  Corinaldesi.  Figure  2  shows  schematically  the  calcula- 

(3) 

tion  of  the  spectral  absorptivities  which  are  carried  out  in  appen¬ 
dices  C  and  D  and  how  they  relate  to  the  other  programs. 

The  population  of  the  initial  level  is  one  of  the  most  impor¬ 
tant  determinants  of  the  strength  of  an  absorption  as  the  population 
in  the  initial  level  can  vary  in  orders  of  magnitude  as  the  temperature 
changes  a  few  thousand  degrees.  If  the  line  in  bound-bound  transition 
arises  from  a  ground  or  very  low  lying  state,  the  population  in  the 
level  is  very  hi^h  in  the  cooler  gas  near  the  walls,  thereby  leading 
to  a  high  spectral  absorptivity.  Examples  for  these  ground  state 
transition  are  the  resonances  lines  of  xenon  lying  in  the  vacuum  ultra¬ 
violet  and  the  yellow  doublet  of  sodium.  The  high  spectral  absorptivity 
will  be  manifested  in  the  arc  discharge  by  ve-"  broadened  and  self 
reversed  emission  bands. 
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Fig.  2-The  information  flow  for  the  calculation  procedures  in  the  appendixes 
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The  shape  of  the  line  determines  how  the  absorption  varies 
witn  wavelength.  We  have  considered  the  lines  in  the  dense  highly 
radiative  plasmas  treated  herein  to  have  a  Lorent-iaa  shape.  The  line 
width,  w,  and  the  shift,  d,  arises  in  the  cases  being  considered  from 
electron  broadening  in  the  center  of  the  arc  (Stark  broadening  theory) , 
and  neutral  particle  broadening  near  the  edge  (Van  de  Waals  or 
resonance  broadening) . 

The  various  broadening  mechanisms  are  calculatable  (in  theory 

at  least)  for  atomic  systems  such  as  xenon.  Figure  3  shows  the 

spectral  absorptivities  for  xenon  for  a  range  of  temperatures  at  the 

same  pressure;  5  atm.  The  program  for  this  is  in  Appendix  C.  The 

spectral  absorptivities  shown  in  Figure  3  include  the  broadened  lines 

and  the  continuum.  The  lines  were  calculated  using  an  intermediate 

(4) 

coupling  program  described  in  the  last  report  and  in  Appendix  B  which 
was  the  extension  to  xenon  of  the  calculation  of  Garstang  and  Van 
BlerkomF^The  continuum  was  that  calculated  by  Schluter ^^Ind  is 
described  in  an  earlier  report  and  in  Appendix  c.  The  upper  limit  on 
the  spectra)  absorptivities  that  may  be  calculated  by  this  program  is 
limited  to  below  those  temperatures  for  which  Xenon  II  becomes 
appreciable  (~  <  14,000°K)  (We  are  modifying  the  programs  to  include 
Xell,  but  this  work  was  not  available  for  this  report.)  The  effect  of 
the  Xell  will  be  to  increase  the  continuum  and  the  emission  in  the  blue 
and  ultraviolet  as  Xell  does  not  have  the  strong  infrared  lines  that 
Xel  exhibits. 

2 . 3  The  Arc  Models 

Two  approaches  were  used  in  the  creation  of  models  for  the 
arc.  In  the  simplest,  the  arc  was  assumed  tc  be  homogeneous  in  tempera¬ 
ture  and  thermal  conduction  to  the  walls  was  neglected.  In  the  more 
complete  approach,  the  absolute  tr-nerature  profile  is  calculated  for 
a  given  electric  field  and  tube  diameter.  From  the  temperature  profile, 
the  spectral  radiance  through  the  arc,  the  radiation  external  to  the 
arc,  and  the  thermal  conduction  to  the  walls  can  be  calculated. 
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The  Inclusion  of  the  free-free  transitions  and  the  more  com¬ 
plete  calculations  of  the  bound-free  transitions  required  more  extended 
application  of  the  quantum  detect  methods  of  Seaton ^ ,  and  Burgess  and 
Seaton^^  .  Peach  (18,19)^  Schluter^^,  and  a  series  of  papers  by 

Biberman  and  Norman  and  their  collaborators  have  indicated  improve- 

(21)  (22) 

ments  using  the  Coulomb  field;  Stewart  and  Rotenberg  arid  McGuire 

have  extended  the  theory  to  the  non-Coulomb  potential  found  in  the  inner 

part  of  the  atom.  The  rise  of  these  methods  required  programs  for  cal- 

(23) 

culating  Whittaker  functions  of  real  and  imaginary  argument.  Appen¬ 
dix  D  presents  ALGOL  translation  of  and  corrections  to  the  FORTRAN  pro- 
(24) 

grains  of  McGuire  .  These  programs  are  to  be  used  as  entry  points  for 
the  numerical  integrations  of  the  radial  matrix  elements  for  bound-free 
and  free-free  transitions. 

The  homogeneous  temperature  models  considered  the  arc  to  be 
a  cylinder  of  uniform  temperature  with  the  arc  diameter  to  be  an  input 
parameter.  The  homogeneous  temperature  models  for  the  pulsed  arcs  were 
found  to  give  reasonable  agreement  with  experiment  if  the  arc  diameter 
were  properly  chosen.  In  the  previous  report,  the  electrical  conduc¬ 
tivity  determined  from  a  .41  cm  bore  tube  was  used  to  calculate  the 
arc  core  diameters  for  the  larger  tubes  at  the  same  temperatures  using 
temperatures  measured  from  a  region  of  the  spectrum  that  was  assumed  to 
be  radiating  like  a  black  body.  The  radiant  emittances  (neglecting  lines) 
calculated  for  the  respective  core  diameters  at  the  temperatures  studied 
were  found  to  account  for  the  major  portion  of  energy  delivered  to  the 
flash  tube. 

The  spectral  radiance  along  a  .8  cm  diameter  for  a  homogeneous 
temperature  discharge  of  12,000°A  and  5  atm.  is  shown  in  Figure  4. 

The  spectral  radiant  emittances  fcr  two  homogeneous  temperature 
cylindrical  arcs  of  .8  and  1.4  cm  diameter  are  shown  in  Figure  5.  The 
radiant  emittance  (i.e.  radiative  flux)  is  shown  in  the  figure.  The 
radiated  power  per  unit  length  is  simply  the  product  of  the  radiant 
emittance  and  the  surface  area  of  the  arc.  The  program  for  the  calcu¬ 
lation  of  the  radiant  emittance  is  part  of  Appendix  C. 
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The  ideal  model  for  the  arc  is  one  in  which  the  characteristics 
of  the  arc  discharge  are  determined  by  the  input  parameters:  the 
electric  field  across  the  arc,  the  geometry  of  the  a^c  container  (i.e., 
the  inside  di-meter  of  the  cylinder),  and  the  pressure  and  composition 
of  the  fill  gas.  But  this  model  requires  the  solution  to  an  integral- 
differential  equation,  called  the  energy  balance  equation,  shown  in  part 
in  Figure  1.  The  result  of  the  application  of  the  equation  is  simply  the 
temperature  distribution  within  the  arc  column. 

We  have  explored  a  few  methods  for  solving  the  energy  balance 
equation  for  non-gray,  nonhomogeneous  plasma.  The  problem  may  be  divided 
into  two  phases:  determining  the  radiative  flux  (i.e.,  radiant  emittance) 
throughout  the  arc  column,  and  using  the  radiative  flux  together  with  the 
transport  properties  to  find  the  temperature  profile  that  satisfies  the 
integral  differential  equation. 

In  this  work,  on  Contract  Nonr  4647(00),  the  radiative  flux  wa_. 

calculated  using  an  exact  integration  (described  in  Reference  2  and  25) 

/on 

and  an  approximate  integration  based  upon  a  technique  of  Sampson'1  1 . 

The  latter  method  was  approximately  ten  times  faster  to  compute  than  the 
former,  and  has  given  good  agreement  for  the  cases  tested  though  it 
does  not  as  yet  include  lines  nor  very  steep  ionization  edges. 

The  temperature  profile  is  calculated  by  using  a  relaxation 
method  applied  to  the  energy  balance  equation.  In  this,  the  radiative 
flux  for  an  assumed  (best  guess)  temperature  profile  is  calculated.  Using 
this  radiative  flux  in  a  time  dependent  energy  balance  equation,  the  con¬ 
tributions  of  the  other  components  of  the  energy  balance  equation  are 
then  calculated.  By  iterating  the  process  over  and  over,  a  temperature 
profile  that  satisfies  the  equation  could  be  calculated  (i.e.,  where  the 
time  dependent  part  equals  zero).  This  is  described  in  Appendix  E. 

A  program  has  been  developed  which  combines  the  exact  and 
approximate  calculations  to  obtain  a  profile.  The  results  of  a  sample 
calculation  together  with  the  program  are  given  in  Appendix  E  also. 

To  provide  further  insight  into  this  calculation,  an  analysis 
procedure  devised  in  other  wor^was  applied  to  the  xenon  arc.  The  complete 
analysis  together  with  calculated  fluxes  and  temperature  profiles  is  in 
Appendix  F. 
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In  all  of  these  calculations,  the  pressure  throughout  the 

discharge  is  considered  to  be  constant;  that  is,  radiation  and  magnetic 

(1  2) 

effects  upon  the  pressure  were  neglected  *  .  The  gas  pressure  before 

the  initiation  of  the  discharge  is  known  (typically,  it  is  150  to  ''00 
torr  for  the  xenon-filled  flash  tubes).  Recent  calculations  have  indi¬ 
cated  that  the  final  pressure  in  the  fully  developed  arc  is  strongly 
dependent  upon  the  temperature  profile.  A  simple  example  will  illustrate 
this  dependence.  Consider  a  simple  two  temperature  confined  arc  in  which 
the  arc  core  is  at  10,000°K,  the  volume  near  the  walls  is  at  1000°K, 
and  neglect  the  electron  contribution  to  the  pressure.  It  is  simply 
shown  that  a  cool  arc  \  olume  of  only  10%  of  the  central  core  volume  will 
reduce  the  arc  pressure  by  a  factor  of  two  from  the  pressure  that  would 

occur  if  the  central  core  occupied  the  whole  volume  as  we  had  assumed 

(1  2  3) 

in  the  earlier  homogeneous  temperature  models  *  *  .  This  consideration 

in  the  pressure  becomes  very  important  for  the  larger  diameters  with 

(4) 

larger  cool  portions 

III.  COUPLING  THE  PUMP  INTO  THE  LASER  ROD 

The  ideal  pump  for  optical  lasers  is  one  in  which  a  major 

fraction  of  the  energy  radiated  can  be  coupled  or  at  least  incident  upon 

the  laser  rod.  There  are  two  general  methods  to  do  this!  che  diffusely 

reflecting  cavity  and  the  focussing  cavity.  We  have  investigated  a 

particular  diffusely  reflecting  cavity,  the  coaxial  laser  pump,  which  is 

well  suited  for  use  with  large  high  energy  pulsed  laser,  and  a  high  image 

quality  focussing  system,  the  spherical  cavity,  which  has  demonstrated 

(28) 

high  efficiency  in  continuously  operating  laser  systems'1 

Figure  6  shows  a  cross  sectional  drawing  transverse  to  the  red 

of  a  current  design  coaxial  lamp.  The  radiation  from  the  plasma  surrounds  the 

laser  rod.  The  MgO  layer  diffusely  reflects  the  radiation  back  toward 

(29) 

the  laser  rod.  As  Whit  .le  and  Skinner  point  out,  in  a  recent  paper, 
the  efficiency  of  the  system  may  be  obtained  roughly  by  comparing  the 
product  of  the  area  and  energy  absorbed  by  the  rod  with  the  sum  of  the 
products  of  areas  and  absorptions  for  the  whole  system.  In  this  analysis 


14 


Dwg.  746A824 


/  Diffusely  Reflecting  Boundary 
/  ic  ter  wall  of  arc  channel) 

—  Arc  Channel 

where  Ky  *  kJV  »  ku(RJ 

-  Arc  Channel  Inner  Wall 


Liquid  with  an  Index  Nearly 
Matching  Arc  Channel 

L  Laser  Rod  in  Which  <wal!  &  ,aser  rod  (n  - 1-5)) 
is  Independent  of  R 


Fig  6— Simplified  representation  of  laser  rod-coaxial  laser  pump  system 
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Figure  7  is  of  a  possible  evolution  of  the  coaxial  1 
which  minimizes  thermal  losses,  and  allows  very  biz,  "“P 
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fr-rtm  e-k  ,  s  ne  cner®al  losses  arisen 
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flecting  surfaces  and  of  thermally  conducting  walls 
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tOCal  »£  ^  ««•  *«  •  »h:VAC  rod  with  external  el «  3 

“  -  -  High  efficiency  system  to  larger  laser  rodTlT 

The  spherical  cavity  has  many  advantages  for  laser  pumping  a-d 
.major  disadvantage,  this  being  the  size  sphere  required  for  a  rod  of  a 
van  ength  (the  diameter  of  the  sphere  should  be  considerably  larger 
n  he  ro  length,.  Some  of  th.  advantages  intlude  a  high  c  an  1 

relatively  remote  reflecting  surface  (less  subject  to  radiation  damage, 
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The  sphere  may  be  fabricated  in  large  sizes  readily  to  optical  tolerances. 
The  rtruli-s  obtained  with  the  spherical  cavity  indicate  promise  for  use 
in  hiph  energy  laser  pumping. 

IV.  RECENT  EXPERIMENTAL  MEASUREMENTS  ON  THE  ARC 

The  experimental  measurements  taken  throughout  this  work  have 
been  directed  towards  improving  and  verifying  the  information  used  in  th 
models.  Early  in  the  work,  we  developed  a  technique  to  measure  the  tem¬ 
perature  of  the  plasma  (More  correctly  we  developed  a  technique  to  measure 
the  lower  bound  of  the  temperature  -  we  will  treat  this  more  thoroughly 
in  the  discussion  in  double  pulsed  arcs).  To  calculate  the  spectral 
properties,  we  needed  to  know  the  pressure  in  the  discharge.  As  we  had 
the  initial  pressure  of  the  gas  (typically  150  torr)  and  the  arc  filled 
the  tube  to  a  first  aDproximation,  we  initially  assumed  a  homogeneous 
temperature.  The  arc  was  measured  from  the  radial  distribution  of  the 
spectral  radiance  in  the  ultraviolet  to  be  relatively  homogeneous  for  a 

large  fraction  of  the  radius  -  though  subsequent  estimates  from  the 

(4) 

electrical  conductivity  discussed  in  the  last  semi-annual  report 
indicated  this  fraction  varied  with  flash  tube  inside  diameter.  Recent 
analysis  of  the  effect  upon  the  pressure  of  the  temperature  profile 
within  the  tube  indicated  that  measurements  of  the  pressure  or  the 
electron  density  were  required. 

The  experimental  work  in  the  past  few  months  has  been  concerned 

with  measurements  of  the  pressure  using  a  piezoelectric  transducer,  a 

Kistler  Instrument  Co.  603A,  and  of  the  spectral  radiance  in  the  vicinity 

of  some  lines  and  in  various  continuum  regions  for  the  arc  under  power 

3  •'  7) 

densities  ovei  1.5  MW/cm  (using  pulse  preionizatior.  '  ).  ?*easurements 

were  also  made  of  the  spectral  radiance  in  the  ultraviolet  (at  3000..)  for 
the  arc  viewed  tide  in  and  viewed  end  on.  These  measurements  from  the 
side  and  end  would  give  a  path  length  difference  of  ten  or  more  depending 
upon  the  tube. 
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The  pressure  measurements  for  the  12.7  and  19  mm  inside 
diameter  tubes  had  shown  a  pressure  considerably  lower  than  that  cal¬ 
culated  for  a  homogeneous  temperature  distribution.  The  transducer  output 
was  very  oscillatory  reducing  the  accuracy.  These  measuremtnts  indicated 
that  there  was  a  cool  volume  near  the  walls  which  could  change  consid¬ 
erably  during  tie  arc  cycle.  A  laser  interferometer  is  being  set  up  to 
measure  the  electron  density  in  the  center  of  cne  arc,  tut  no  rinal 

results  are  available  for  this  report.  The  interferometer  is  similar  to 

(8) 

that  of  Gerardo  et.  al.  in  which  the  flash  tube  is  inside  a  reference 
cavity  which  is  in  series  with  the  laser  cavity.  Measurements  of  the 
electron  density  with  the  interferometer  together  with  those  of  temper¬ 
atures  will  characterize  more  'ailly  the  arc  plasma. 

The  measurements  of  the  temperature  in  the  double  pulsed  arc 
required  finding  spectral  regions  in  which  the  emission  from  the  arc  was 
not  self-reversed  but  still  optically  thick.  The  8231. 6A  line  of  Xel 

3 

was  suspected  to  be  self-reversed  at  power  densities  above  about  .1  MW/cm  . 

This  was  confirmed  using  the  high  resolution  rapid  scanning  spectrometer^ 

Calculations  of  the  spectral  radiance  (using  the  program  in  Appendix  C) 

for  a  temperature  profile  similar  to  that  expected  in  the  arc  showed  self 

o 

reversal  in  the  immediate  vicinity  of  the  8231.6  A  lines. 

The  spectral  radiances  measured  for  the  lines  and  the  continuum 

3 

are  shown  in  Figures  9  through  12  for  power  densities  over  1.5  MW/cm  . 

In  Figures  9  and  10,  the  spectral  radiance  is  plotted  in  terms  of  the 

equivalent  black  body  temperature.  The  power  is  the  input  power  dissipated 

in  the  lamp  per  unit  volume.  The  flash  tube  diameters  ranged  from  .41  cm 

to  1.94  cm.  The  tubes  were  preionized  with  a  .5  ms  pulse  (100-600  uF  in 

series  with  100  yHy) .  The  main  pulse  was  about  50  us  long  (50-100  uF  and 

the  residual  inductance  ''■3uH) .  The  main  pulse  was  switched  in  by  means 

of  an  ignitron  using  the  circuit  shown  in  Reference  6.  The  self  reversal 
o 

of  the  8231.6  A  line  of  Xel  was  considered  to  be  the  primary  reason  for 

the  reduction  in  the  spectra]  radiances  for  the  larger  diameter  flash 

tubes.  Figure  13  gives  the  spectral  radiance  distribution  for  an  EG  &  G 

Inc.  FX-52  flash  tjbe  under  single  and  double  pulsed  conditions.  The  data 

o 

was  taken  with  a  spectral  width  of  about  1  A  at  the  wavelengths  indicated. 
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Peak  Power  Density,  kW/cm 

Fig.  10- Black  body  U  nperature  corresponding  to  same  spectral  radiance  as  13  mm  flash  tube 

as  a  function  of  input  power  density 
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Fig.  13  Spectral  radiancedistribution 


The  measurements  of  tha  spectral  radiance  at  3000  X  on  flash 
tubes  viewed  side  on  and  viewed  end  on  yielded  the  theoretical  ratio  of 
the  length  of  the  tube  to  diameter  of  the  tube  for  lower  energies  only 
when  the  tungsten  electrodes  were  moved  back  out  of  the  arc.  The  ratio 
again  deviated  at  the  higher  currents  as  shown  in  Figure  14.  This  may 
arise  in  part  from  further  tungsten  injection,  but  also  may  be  due  to  the 
plasma  becoming  thick  for  the  flash  tube  length  being  viewed. 

We  have  not  yet  made  a  detailed  comparison  of  the  calculated 
radiances  with  those  measured  experimentally.  This  requires  a  more 
accurate  knowledge  of  the  pressure  and/or  the  electron  density  than  we 
currently  have.  We  shall  be  making  these  comparisons  after  the  laser 
interferometer  is  completed.  To  compare  the  spectral  radiances  observed 

in  Figure  13  with  the  calculated  values  in  Figure  4,  multiply  the  values 

8  2 
in  Figure  13  by  10  to  obtain  the  values  in  Figure  4  (W/cm  Ster  cm) . 

The  single  pulsed  values  correspond  with  those  for  the  plot  in  Figure  4; 

the  exact  agreement  requiring  better  values  of  pressure  than  we 

currently  have. 
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Fig.  14— The  spectral  radiance  viewed  at  3000  A  end  on  and  side  on 
for  a  13  mm  diameter  arc  tube 


V.  SUMMARY  AND  CONCLUSIONS 

This  report  presents  a  description  of  the  work  carried  out  ca 
Contract  t >nr  4647(00)  towards  the  creation  of  quantitative  models  for 
the  pulsed  arc  discharge  that  is  used  for  the  optical  pumping  of  high 
energy  solid  state  lasers.  The  primary  results  of  these  studies  are 
computer  programs  for  the  ~alculations  of  the  physical  properties  of  the 
arc  plasma,  using  xenon  as  the  working  gas,  and  for  determining  the 
balance  within  the  arc  between  radiation  external  to  the  arc  and  thermal 
losses  to  the  walls.  The  experimental  results  obtained  since  the  last 
semiannual  reports  are  also  described.  The  experimental  studies  have 
been  concerned  with  measurements  of  the  pressure  within  the  arc,  and  the 
characterization  of  the  arcs  at  very  high  Dower  densities. 

The  models  for  the  arc  developed  in  this  work  are  at  present 
only  semi-quantitative.  They  possess  the  main  features  of  the  arc, 
exhibiting  the  saturation  of  the  radiation  in  the  infrared,  the  strong 
dependence  upon  power  density  of  the  radiation  emitted  in  the  ultraviolet; 
and  the  nearly  homogeneous  temperature  distribution  that  have  been 
observed  in  these  arcs.  To  improve  the  accuracy  of  the  models,  more 
complete  and  accurate  measurements  of  the  physical  prope-ties  need  to  be 
made.  Through  the  use  of  the  model  calculations,  quantitative  compari^ono 
can  be  made  between  the  power  input,  the  average  electrical  conductivity, 
the  spectral  absorptivities  in  various  regions  of  tV  nontinuum  and  in  the 
vicinity  of  selected  lines,  and  of  the  pressure  and  electron  density  to 
make  a  consistent  picture. 

The  full  implications  of  these  models  are  still  being  explored. 
The  models  are  allowing  a  ■  ew  insight  into  the  behavior  of  the  big  .  density 
arc.  The  arc  models  help  to  remove  the  high  density  plasma  -  the  light 
emitting  ~lasmas  -from  the  "art"  area  of  technology. 

We  hope  that  che  work  described  in  this  series  of  reports  will 
assist  the  designer  of  optically  pumped  lasers  in  the  development  and 
improvement  of  this  area  of  laser  technology. 

Through  this  work,  we  hope  that  v.e  are  bring  to  bear  on  an  old 
problem—  light  sources — some  of  the  new  techniques  in  both  theory  and 
experiment  that  were  developed  by  a  wide  variety  of  people  working  in 
areas  far  removed  from  the  optical  pumping  of  laserr 
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APPENDIX  A 


A  PROGRAM  FOR 
OF 


THE  CALCULATION  OF  THE  ELECTRICAL  AND  THERMAL  CONDUCTIVITY 
XENON  AS  A  FUNCTION  OF  TEMPERATURE  AND  PRESSURE 


by 


R.  S.  deVoto 
Stanford  University 


This  program,  ELEC/THERM,  calculates  the  various  components  of 
the  electrical  and  the  thermal  conductivities  of  xenon  as  a  function  of 
pressure  and  temperature  under  local  thermal  equilibrium  conditions. 

The  program  will  calculate  these  quantities  for  argon  also. 

The  input  data  for  argon  is: 

1.  "ARGON",  number  of  species  -  1  (we  used  2) 

2.  initial  temperature,  temperature  step  size,  maximum  temperature 
desired 

3.  number  of  pressures;  the  pressures  in  atmospheres 

4.  PHIOA  (3.23  @  4),  RHOA  (.224),  PHIOI  (5.38  @  6),  RHOI  (.1965), 
AC  (25.615),  BC  (1.1960),  HIER  (must  be  Oj . 

5.  electron  mass  (9.1091  at  -28),  molec.  wt.  of  atom  (39.944), 
molec.  wt.  of  first  ion  (39.944),  molec.  wt.  of  any  other  ions 
used. 

For  Aenon  the  data  required  is  similar  to  the  ARGON  program. 

1.  "XENON",  2, 

2.  initial  temperature,  step  size,  maximum  temperature  desired 

3.  number  of  pressures,  the  pressures  (in  atmospheres) 

4.  3.11  @  6,  .208,0,0,  25.639,  .99752,0,  9.109  @  -28,  131.3,131.3. 
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APPENDIX  B 


CALCULATION  OF  BOUND-BOUND  TRANSITION  PROBABILITIES 

by  D.  Sun 

The  overall  purpose  of  the  system  or  programs  being  developed 
is  to  represent  the  spectrum  of  an  atomic  species  in  a  plasma,  including 
all  lines  and  portions  of  the  continuum  which  carry  a  significant  portion 
of  the  energy.  Hence  it  will  oe  necessary  to  calculate  large  numbers  of 
transition  probabilities  for  bound-bound,  bound-free,  and  free-free 
transitions. 

Calculation  of  the  bound-bovnd  case  to  be  discussed  here  is 
fairly  complete;  the  other  two  types  of  transitions  have  not  been  ad>’ed 
yet.  The  program  to  be  described,  presently  called  GIANT/FIASCO,  con¬ 
sists  of  four  sub-programs,  (following  the  format  of  E.  Corinaldesi) 
between  which  information  is  transferred  by  three  disc  files,  READFl, 
RE.J)F2,  and  READF3.  Programs  1  and  2  are  run  together,  and  accepting 
as  their  input  the  quantum  numbers  and  energies  of  the  various  energy 
levels  of  a  species  as  they  are  given  in  Moore's  tables.^  Program  2 
prints  out  transition  probabilities  for  all  transitions  allowed  by  the 
select!  :  rules  within  the  range  considered.  These  are  also  placed  on 
the  disc  READF2  where  they  may  be  tapped  by  Program  ?.  or  by  an  auxiliary 
program  FIASCG/RESULTS ,  which  merely  prints  out  the  results  In  a  number 
of  convenient  forms.  Program  3  calculates  the  Stark  shift  and  broaden!,:;; 
of  each  level,  while  Program  4  combines  these  into  the  shift  and  broaden¬ 
ing  of  each  spectra*,  line.  However,  this  paper  is  concerned  primarily 
with  Programsl  and  2,  the  calculation  of  bound-bound  transition  probabili¬ 
ties  . 

2 

The  heart  of  the  matter  is  the  calculation  of  1'ijs  Izlii'  | 

n 

where  ,  and  s{i  are  the  two  states  between  which  a  transition  may  occur 
n  m 

Moore,  Atomic  Energy  Levels’,  National  Bureau  of  Standards,  Circular 
467,  U.S.  Government  Printing  Ifrirc,  Washington,  D.C. 

n 

‘F.  Corinaldesi,  to  be  published. 
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and  zmr  cos  0  is  the  electric  dipole  operator.  (Quadrupele  and  magnetic 
dipole  transitions  are  ignored.)  Although  tj>n,  v  ,  and  z  actually  involve 
all  the  electrons  in  the  atom,  in  practice,  it  is  possible  to  use  the 

one-electron  wave  functions  ct  tht  cter  electron  which  makes  ti.e  transi- 

,  •;  3 

tion  (Condon  and  Shortley^  Section  6  ).  Now  each  energy  level  consists 
cf  2J+1  states  with  different  z-components  of  angular  momentum.  For 
transitions  from  n  to  m,  it  is  necessary  to  sum  the  square  of  the  matrix 
element  for  all  of  the  m-states,  and  average  it  over  all  of  the  initial 
n-states.  The  oscillator  strength  f  is  J-  times  this  quantity,  and 
the  following  relations  hold:  n 

(2  J  +  l)f  -  gf  -  (2  J  +  l)f 
n  nm  m  mn 

-o  _  _1 

A  *  0.667  x  v  where  v  is  in  cm 
nm 


(2  J  + 
n 


1)A. 


nm 


(2  + 

m 


1)A 


mn 


The  solutions  to  Schroedinger ' s  equation  separate  into  two 
parts,  giving  the  form  Y(0,  <p) .  R(r)  is  known  as  the  radial  wave 

function  while  Y(0,  <D  is  the  angular  part  Hence  the  problem  separates 

JUS 

R^^rdr,  and  an  averaging  and  summing  process 
which  is  carried  out  using  Racah  coefficients  of  the  quantum  numbers  cf 
the  initial  and  final  states. 


THE  ANGULAR  PART  (COUPLING  SCHEMES ) 


When  an  atom  with  more  than  one  electron  is  being  considered, 
there  are  several  ways  to  add  the  multiple  orbital  and  spin  angular 
momentums.  For  example,  one  might  add  all  the  orbital  angular  momentum 


3 

E.  Condon  and  G.  Shortley,  Theory  of  Atomic  Spectra,  Cambridge  University 
Fress,  London  (1935)  hereafter  reft  rred  toasCondon  and  Shortley  and  the 
Section. 
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operators  to  get  a  total  orbital  angular  momentum  operator,  and  do  the 
same  with  the  spin  operators.  Then  the  total  orbital  mgular  ..■omentum 
operator  and  the  total  spin  angular  momentum  operator  might  be  added  to 
give  a  total  agular  momentum.  If  the  actual  states  of  the  atom  closely 
approximate  the  eigenstates  of  the  total  orbital  and  spin  angular  momentum 
operators,  then  the  atom  is  said  to  be  LS  coupled,  the  L  and  S  repre¬ 
senting  good  quantum  numbers  of  the  operators.  One  way  to  express  the 
states  is  to  put  the  quantum  numbers  of  two  operators  in  brackets  to  the 
left  of  the  quantum  number  for  the  sum  of  the  two  operators.  An  LS  state 

is  written  [S,  (L  ,  H)L)J>,  where  L  is  the  total  L  of  the  inner 

1  core  core 

electrons,  i  is  the  orbital  quantun  number  of  the  outer  electron,  and  .1 
is  the  total  angular  momentum.  A  form  of  coupling  commonl'  found  in  the 
rare  gases  is  coupling,  given  by  U(Jcore*  OK,  S]J>,  where  Jcore  is 
the  total  angular  momentum  of  the  core,  s  is  the  spin  of  the  outer 
electron,  and  K  is  the  intermediate  quantum  number. 

The  average  and  sum  over  the  initial  and  final  states  for 

2 

|<i|>  IzU  >1  is  easily  seen  to  be  equ.i  to  one  third  of  the  same  average 
n  m  „ 

and  sum  for  |  <>£  |r|'<l/  >1  .  r  turns  out  to  be  a  class  T  operator  with 

n  3 

respect  to  J  (See  Condon  and  Shortley,  8  ).  Here  it  is  convenient  to 

introduce  a  quantity  <j;Tij'>  which  can  be  related  to  the  components 

of  the  matrix  element  jm  j  T |  j 'm' >  but  is  itself  independent  if  m  and  m' , 

3 

the  z-components .  (See  Condon  and  Shortley,  9  ).  Now  if  [P,J^]  »  0 
and  P  is  class  T  with  to  J  =  then  the  ratio 

i<jlj2iF:"jlj2,>|2 

3 

is  given  by  Equation  11  8  of  Condon  and  Shortley,  or  by 
(2j  '  +  1) J 2> (2j 2  +  l)(2j2’  +  L)W2(j2:j,j2,,j’;j1)/  T  (j,j’)  where 
j2>  =  max(J2,j2'),  W(j2,j,j2',j,;j1>  is  the  Racah  coefficient,  and 
3  (j,j')  is  defined  by: 
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X  (j,j+l)  -  <J+l)(2j+3) 
X  (j.j)  “  j(j+l) 


X  “  j(2j-l)  See  Condon  ar.d  Shortley,  7^5 

In  the  program  there  are  two  procedures,  XI  and  W2M0D,  with  XI (j  * j ' )  = 
X  and  with  W2M0D(j  2 ,  j  ,  j  2' ,  j  ' ,  /XI  ( j  ,  j  ' )  equal  to  the  above 


ratio. 


Returning  to  LS  and  jl  coupling,  we  have: 


I  <  [S ,  (L  ,  £)L]Jfri  [  S ,  (L  ,£')L']J'>|  = 
'  core  :  -  core  1 


W2M0D  Jj.LJjJ  '  XSJ  ,  i)Ljri(L 

XI(J, J! )  1  lLcore»^L  r:  Uco.e,£  1 


W2M0D(L,J,L'  ,J'  ,S)  W2M0D ,L'  ,L  ) 


core 


XI(J.J')  XI(L,L  ) 


|<£ir  •;  £ ’>  |  2 


for  LS  coupling,  and 

|<[(J  ,£)K,s]J:r j [J  , £ ' )K’ ,s] J ' > | 2  - 

1  core  core  1 


W2M0D (K, J., K  , J — j  £)Kir{(J  £'w»- 

XI(J ,  J ' )  1  Ucore’*;K:rt(Jcore’£  ' 


W2M0D(L,J,K',J\s)  W2M0D(£,X,£.,,K' ,J  )  „ 

I  <  3  ”r  J  C  ’>  |  2 

XI(J,J’)  XI(K,K')  1  ♦  •  1 


for  j£  coupling. 

Summing  over  all  z-components  of  J  for  both  initial  states  and  final 

states  and  dividing  by  (2J+1)  is  equivalent  to  multiplying  by  XI(J,J'). 

4  3 

(See  Condon  and  Shortley,  Equations  7  5  and  9  11) 
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|e  r  t'>|2  =  - — —  [fR  R’  rdr]2 

4  i  -1  0 

> 

(R  and  R'  are  radial  functions  although  the 

quantum  numbers  pertaining  to  the  radial  part  have  not  been  written  in 

2 

the  abo”e  equations.  )  ihe  final  result  for  3 1 <n j  z |m>  j  summed  over  final 
states  and  averaged  ever  initial  states  is: 

W2M0D(L,  J ,  L ' ,  J  ' ,  S)  W2M0D(S,,L,i' ,L' ,L  oj.e)  [/" R  R'  rdr]2 
XI (L, L ' )  (2i+l)(2t'+l) 

for  LS  coupling,  and 

W2M0D(L,  J ,K' , J ' , s )  W2M0D ( S. ,  K ,  £ r  ,K'  ,J  )  [/"rR'  rdr]2 

core  6 

XI(K,K')  (2£+l)(2£'+l) 

When  these  quantities  are  multiplied  by  (2J+1)  and  summed  over 
all  possible  J  and  J'  in  the  initial  and  final  ujultiplets,  the  result  is 
in  the  notation  used  by  *<ohrlich^  and  Griem,^  the  7H  standing  for 
"multiplet  s  cac  )6o»)  is  then  (2J+1)  times  the  above  quantities, 
the^C  standing  for  "line".  er  )  obviously  was  defined  to  satisfy  the 
rule  (<£)  =  1. 

A  sum  rule  derived  from  first  principles  states  that  v|<n|z|ra>| 
summed  over  all  discrete  final  states  and  integrated  over  continuum  final 
states  should  equal  the  Rydberg  constant  in  the  same  units  as  v.  In  the 
program  print  out  the  quantity  S3  or  SUM3  is  the  sum  over  discrete  states 
of  3  v|<n'z]m>|  with  v  in  cm  .  It  should  approach  329,212. 

V  Rohrlich,  Ap.  J.  129,  441,449  (1959). 

"*H.  Griem,  Plasma  Spectroscopy,  McGraw  Hill  Book  Company,  New  York  (1964). 
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Intermediate  Coupling 


fi; 


Most  elements  do  not  fit  a  well-defined  coupling  scheme  as 
nicely  as  might  be  desired.  In  other  words,  their  states  are  mixtures 
of  pure  states.  No  matter  what  basis  is  chosen — LS,  jl,  jj,  etc. — 
the  actual  states  are  linear  combinations  of  the  basis  states,  and  a 
transformation  matrix  between  the  actual  states  and  the  pure  states 
exists. 

For  the  intermediate  coupling  (IC)  program,  an  I.S  basis  was 
used,  and  the  elements  of  the  transformation  matrix  between  states  of 
different  configurations  were  assumed  to  be  zero.  For  each  type  of  con¬ 
figuration  involved,  a  theoretical  matrix  expression  for  the  electro¬ 
static  and  spin-orbit  perturbations  to  the  Hamiltonian  was  generated  in 
terms  of  F^,  GQ,  F2,  G2,  and  C ' •  The  electrostatic  parameters,  FQ, 

G  ,  Fy,  and  G~,  appeared  only  on  the  diagonal  for  the  LS  basis  (See 
Condon  and  Shortley,  1  ),  while  the  two  spin-orbit  parameters  £  and  £ ' 

appeared  in  many  off-diagonal  positions  (Condon  and  Shortley,  l^) . 
However,  no  cross  terms  exist  between  states  of  different  J,  implying 
that  J  is  a  good  quantum  number  in  all  coupling  schemes. 

If  the  theoretical  expression  is  correct  and  if  it  is  valid 
to  ignore  configuration  interaction,  then  it  should  be  possible  to  choose 
values  for  the  parameters  such  that  the  diagonalized  matrix  corresponds 
closely  to  the  observed  energies,  in  spite  of  the  fact  that  the  number 
of  energies  in  a  configuration  (the  number  of  conditions  that  can  be 
imposed  on  the  parameters)  is  larger  than  the  number  of  parameters.  In 
practice,  the  parameters  are  adjusted  to  give  a  least  squares  fit  between 
the  energies  calculated  by  diagonalization  and  the  observed  energies. 

This  is  done  by  using  trace  conditions  from  the  blocks  of  different  J 
within  a  configuration  (There  are  no  interaction  terms  between  states 
of  different  J.)  to  set  several  of  the  parameters  initially.  As  a  first 
guess,  the  less  significant  parameters  are  set  equal  to  zero.  The  matrix 
is  then  diagonalized  and  partial  derivatives  of  the  eigenvalues  with 
respect  to  the  parameters  are  calculated.  A  multiple-dimensional  Newton's 
method  is  applied  to  readjust  the  parameters.  The  process  is  repeated 
until  the  sum  of  the  squares  of  the  energy  differences  changes  by  less 
than  0.1%  from  its  previous  value. 
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The  matrix  elements  for  z  between  various  LS  states  are  easy 
to  calculate,  although  phases  are  somewhat  of  a  problem.  (See  Condon 
and  Shortley,  4  ^ ) .  Using  the  transformation  matrices,  calculation  of 
oscillator  strengths  is  straightforward. 


THE  RADIAL  INTEGRAL 

Schroedinger ' s  equation  for  the  radial  wave  function  is: 

2 

R  =  R  /r  =  P  / r  +  [-2V(r)  -  +  c]R  =  0 

d  iU  3r2  r2 

where  r  is  in  Bohr  radii 
and  e  is  in  Rydbergs 

Both  of  the  methods  used  in  these  programs  for  solving  this 
equation,  the  Bat.es-Damgaard  approximation^  and  the  scaled  Thomas-Fsrmi 
method,'7  are  semi-empirical — that  is,  e  is  set  equal  to  the  observed 
energy  of  a  known  level.  The  difference  between  the  two  approaches  lies 
in  the  potential  V(r)  that  is  chosen. 

The  Bates-Damgaard  or  Coulomb  Approximation 

As  a  first  approximation  one  guesses  that  the  potential  in 
which  the  outer  electron  moves  is  simply  a  Coulomb  field  due  to  a  nucleus 
of  charge  +Z  shielded  by  NEL-1  electrons  very  close  to  the  nucleus. 

Q 

Thus  V(r)  =  -—where  C  =  Z-NEL+1  is  the  effective  core  charge.  (C  =  1 
*  C^ 

for  neutral  atoms.)  Letting  e  =  — -  the  differential  equation  is  now: 

n~ 


3r 


+ 


, 2C  1(1+1) 

1  r  ”  2 

r 


~  ]  R  =  0 
n 


^L.  Bates  and  A.  Damgaard,  Phil.  Trans.  Roy.  Soc .  (London)  242A,  101  (1949). 
J*  Stewart  and  M.  Rotenburg,  Phys.  Rev.  140 .  A1508  (1965). 
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There  are  two  boundary  conditions,  that  R-+0  at  the  origin  and 

infinity.  They  can  both  be  met  only  when  n  is  an  integer.  In  complex 

atoms  n  is  generally  non-integral,  and  one  condition  must  be  dropped. 

Since  the  major  contribution  to  |  R  R  rdr  comes  from  large  r  for  most 

J0  n  m 

transitions,  it  is  reasonable  to  forget  the  condition  at  the  origin. 

The  wave  function  then  blows  up  at  zero  and  is  not  normalizable. 

A  mathematical  trick  is  applied  at  this  point.  An  asymptotic 

cr 

expansion  for  large  r  is  generated  by  making  the  substitution  R=ue  n 
and  finding  a  recursion  relation  for  u. 


ci 

n 


<^f>n  a  + 


a 

y  ] 

til  rl 


a^  is  given  by  a  recursion  relation  in  reference  6. 

As  long  as  the  £  contains  a  finite  number  of  terms,  the  condition  at 
infinity  holds  and  the  one  at  zero  does  \;c But,  a  quick  check  shows 
that  the  differential  equation  has  an  irregular  singularity  at  infinity. 
Therefore,  a  power  series  expansion  in  (1/r)  does  not  necessarily  con¬ 
verge  to  a  solution.  In  fact,  the  sum  diverges  here.  In  spite  of  this, 
if  the  series  is  terminated  at  the  proper  term,  it  gives  very  good  fit 

to  the  actual  solution  for  r>2  or  3  (depending  upon  n  and  1,  of  course). 

/«! 

R  R  rdr  is  now  an  integral  of  an  exponential  term  times  a 
o  n  m 

power  series  in  r,  and  may  be  expressed  as  sum  of  gamma  functions. 
Because  the  function  does  not  have  to  be  evaluated  numerically  at  a  number 
of  points,  very  rapid  computer  calculations  are  possible.  (The  normali¬ 
zation  is  also  given  analytically.^)  To  overcome  the  problems  of  the 
singularity  of  the  origin  and  of  the  diverging  series,  all  powers  of  r 
less  than  2  in  the  quantity  R  R  r  are  removed.  These  terms  do  not  exist 
physically  since  the  actual  potential  is  not  a  Coulomb  field  near  the 
origin.  For  the  details  of  the  theory  and  calculation,  see  References  6 
and  8  • 


i.  Coriraldesj.  and  E.  Ceil,  Computer  Program  for  Pates-Damgaard  Integrals, 
Westinghouse  Research  Laboratories  Scientific  Paper  66-1C1-EPLAS-P4  (196b). 
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The  Thomas-Fermi  Method 


The  Thomas-Fermi  approach  utilizes  a  somewhat  more  realistic 
potential.  Instead  of  assuming  that  the  N  -1  shielding  electrons  are 
clustered  within  an  infintesimal  distance  of  the  nucleus,  it  predicts  a 
distribution  of  the  shielding  electrons  based  on  Fermi-Dirac  statistics. 
The  fact  that  the  potentials  formed  in  the  above  manner  bring  the  wave 
functions  much  closer  to  fitting  the  boundary  condition  at  the  origin 
justifies  the  use  of  the  method.  Moreover,  the  fit  can  be  made  perfect 
by  linearly  expanding  or  contracting  the  potential  in  space  by  a  scaling 
factor  close  to  one.  The  latter  is  the  scaled  Thomas-Fermi  method. 

The  easiest  way  to  understand  the  method  is  to  look  at  the  net 
charge  contained  within  a  sphere  of  radius  r  as  a  function  of  r. 


The  contained  charge  must  obviously  go  to  Z  as  r-K),  while  for  r>ro — 

that  is,  for  r  outside  of  the  core — the  contained  charge  is  exactly  C, 

as  in  the  Bates-Damgaard  method. 

( 


containing  I  Z  [  4>  (— )  +  —  ] 
•  ~  \  u  r 

charge  ' 


V(r)  = 


-Z[*(“)  +~] 

J3L 


!  -  C/r 


r<r 


r>r 


where  <Ko)  *  1  ,  4>  (x  )  =  0  ,  $ '  (x  )  - - — 

O  O  X  L 

O 


with 


_ o 

o  u 


According  to  Fermi-Dirac  theory: 


u  =  0.8853  Z 


."1/3 


3/2  -1/2 

~  j  =  4>  x 

dx' 
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At  first  glance  there  appears  to  be  one  extra  boundary  condition 

on  the  differential  equation  for  $.  However,  xq  has  not  been  specified, 

and  all  three  conditions  are  needed  to  uniquely  determine  <£ .  The  equation 

must  be  solved  numerically,  but  without  x^  the  calculations  are  quite 

extensive.  Luckily,  the  work  has  already  been  done,  and  a  polynominal 

fit  for  the  value  of  x  as  a  function  of  (C/Z)  has  been  made.^  The 

o 

present  J1  and  IC  programs  merely  solve  for  x  using  the  polynomial  and 

°  r 

then  generate  $  by  starting  at  xq  with  $(xq)  *  0  and  $'(xq)  - - —  .  <p 

automatically  becomes  one  at  zero.  This  calculation  need  only  be  performed 
once  for  each  atomic  species,  for  $  is  completely  determined  by  C  and  Z. 

's  equation  must  also  be  done 
s  possible  because  V(r),  which 

S 

y.  Two  transformations  are  made:" 
R(r)  *  R(y2) 
dr  ■  2ydy 


[ 2V (y 2 )  -  cj  q 

J 

This  is  solved  by  Numerov's  method  which  has  an  error  in  each  step  on  the 

6  1/2 

order  of  h  (h  •  step  size).  In  the  program  h  was  yQ/32  where  yQ  =  rQ 

The  integration  was  started  at  ay  quite  a  bit  larger  than  y^ 

by  evaluating  the  asymtotic  expansion  of  the  Batts  and  Damgaard  (BD) 

function  at  two  points  separated  by  h.  Care  must  be  taken  that  y  is 

large  enough  for  the  asymtotic  expansion  to  be  sufficiently  accurate  and 

for  the  tail  of  the  wave  function  to  be  negligible  beyond  y.  Still,  the 

starting  point  must  not  be  so  far  out  that  Numerov’s  method  drifts  from 

the  actual  solution.  The  first  term  of  the  asymtotic  expansion  is  of 

— y  2 

the  form  e  ,  and  such  functions  do  exhibit  considerable  drift  when 
Numerov's  method  is  started  too  far  out.  (The  higher  derivatives  of 


e  solution  to  Schroedinger 
numeric  No  analytical  solution  i< 
depends  c  <p,  is  on^y  known  numericall; 


R(r)  •  p(y2)  ,  r  *  y2  * 


and  p(y)  -  y1^2  q(y)  , 


giving: 


ay 


i(  Z+  1)  4  3 


*y 


i  +  Ay 
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e  'r  are  very  large  compared  to  the  function  itself  when  y  is  large.) 
in  this  program  the  integration  was  started  at  the  point  where  the  first 
term  of  the  Bates  and  Damgaard  expansion  has  dropped  to  about  1/3000  of 
its  maximum  value. 

If  the  scaled  Thomas-Fermi  method  is  being  done,  the  integration 

is  stopped  at  one  step  from  the  origin  and  the  number  ~>f  nodes  in  the 

wave  function  is  checked  to  see  if  it  is  equal  to  n,  -  t  -  1.  n. 

int  mt 

is  the  given  principal  quantum  number,  not  the  effective  n  calculated 

g2 

from  e  =  — j  .  n.  is  always  an  integer.  If  the  number  of  nodes  is 
n 

not  correct,  then  the  direction  in  which  the  scaling  [actor  must  be 
changed  is  immediately  determined.  On  the  other  hand,  if  it  is  correct, 
then  the  ratio  of  the  two  pcints  closest  to  the  origin  is  checked  against 
the  ratio  predicted  by  a  special  asymtotic  expansion  for  small  y,  obtained 
from  a  recursion  relation  for  the  approximate  equation: 


2  .  2 
3y 


From  this  the  direction  cf  the  next  change  is  found.  Tne  scaling  factor 

steps  in  the  proper  direction  by  0.05,  to  either  1.05  or  0.95,  and  the 

diff.  eqn.is  resolved  with  the  appropriately  changed  potential^.  The 

scaling  factor  continues  to  move  in  one  direction  until  the  value  which 

gives  a  boundary  fit  is  overstepped.  Then  the  step  size  is  cut  in  half 

each  time  a  new  step  is  made.  In  this  manner,  the  interval  in  which  the 

correct  value  falls  is  bisected  with  each  step.  The  iteration  goes  on 

until  the  step  size  is  less  than  1/2000  or  until  the  ratio  of  the  two 

points  next  to  the  origin  is  within  one  percent  of  its  predicted  value. 

In  the  rare  gases  the  gr.  nd  state  has  too  small  an  n^.^  to  be 

put  into  the  Bates  and  Damgaard  (BD)  method.  Instead  the  wave  function  of 

the  ground  state  is  calculated  numerically  by  the  scaled  Thomas-Fermi 

methoa.  Therefore,  when  the  ground  state  is  involved  the  /  R  R  rdr  is 

o  n  m 

done  numerically,  implying  that  the  wave  functions  of  ail  upper  states 
that  make  transitions  to  the  ground  state  must  be  known  numerical lv . 

One  way  of  getting  the  upper  state  wave  function  would  be  to  sum  BD 
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asymtotic  expansion  at  each  point  desired.  However,  it  is  quicker  to 
simply  solve  Schroedinger's  equation  numerically  since  Numerov's  method 
requires  only  a  small  number  of  operations  at  each  point,  whereas  summing 
the  asymtotic  expansion  is  somewhat  involved.  As  long  as  the  differential 
equation  is  being  solved  numerically,  it  is  just  as  easy-and  certainly 
more  accurate— to  use  the  Thomas-Fermi  potential  in  place  of  a  pure 
Coulomb  potential.  This  constitutes  the  unsealed  Thomas-Fermi  method 
and  has  roughly  a  factor  of  ten  advantage  in  computer  time  over  the 
scaled  TF  method.  In  general,  wave  functions  generated  this  way  will 
blow  up  at  the  origin  and  must  be  chopped  off  near  r  -  0.  The  program 
imposes  the  form  r  \  the  first  term  of  the  expansion  of  the  wave 
function  for  small  r,  when  the  centrifugal  potential  term  16  £ii±.U  -13 
begins  to  dominate.  4y^ 

GIANT/FIASCO  first  solves  for  <f>  and  then  generates  a  scaled 
Thomas-Fermi  wave  function  for  the  ground  state  and  unsealed  Thomas-Fermi 
wave  functions  for  all  other  states  up  to  the  point  where  n,ff  is  greater 
than  the  ngff  of  the  ground  sfate  by  more  than  DNEFF  an  arbitrary  limit 
fed  into  the  program.  (Transitions  with  a  difference  in  n  greater 
than  DNEFF  are  not  considered.  This  is  true  of  non-gromu^state  transi¬ 
tions  as  well.)  All  generation  of  any  kind  of  TF  wave  functions  is 
performed  by  an  internal  procedure  called  THOMASFERMI.  The  integration 
°f  f”RnRm  rdr  is  carr*ed  out  numerically  when  both  wavs  functions  are 
known  numerically  and  is  done  by  the  Bates  and  Damgaard  method  in  other 
cases.  Of  course,  it  is  not  done  at  all  if  the  difference  in  the  n  's 
is  larger  than  DNEFF.  ef* 
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There  are  several  arbitrary  conditions  dealing  with  the  Thomas- 
Fenni  method  specifically  that  are  built  into  the  program  and  may  require 
modification  in  the  futu'e.  They  are: 

1)  the  step  size  h  for  solving  differential  equatior,s--it  is 
x  xq/4096  for  4>,  and  yQ/32  for  the  wave  functions 

2)  the  conditions  lor  stopping  the  iteration  in  the  scaled 
TF  method 

3)  the  method  of  chopping  off  the  singularity  at  the  origin 
in  the  unsealed  TF  process 

4)  choosing  the  y  outside  of  yQ  at  which  to  start  solving 
Schroedinger 's  equation 

5)  deciding  which  transitions  should  be  done  by  the  BD  method 
and  which  should  be  done  numerically 

6)  deciding  which  levels  should  be  scaled  TF  and  which  should 
be  unsealed  TF. 
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SYMBOL  TABLE 


Z:  the  atomic  number 

NEL:  the  number  of  ele -.crons 

^  th*  COre  char8e  seen  ‘y  'he  outer  electron  C-2-NEL+1 

the  Rydberg  constant  in  reciprocal  centimeters  R., 109737.  31 

an  index  running  over  the  energy  levels 
orbital  angular  momentum  of  the  outer  electron 


R: 

N: 

L: 

J: 

T: 


total  angular  momentum  of  the  atom 
energy  of  the  level  in  cm' i 


TINF:  energy  limit  of  the  level  in  cm-1 
In  the  LS  program  only 

MOLT:  multiplicity  -1ST1  where  S  is  the  total  spin  of  the  atom 
total  orbitai.  angular  momentum  of  the  atom 
P-*  parity  of  the  state 
In  the  jl  and  IC  programs  only 

CORE:  total  angular  momentum  of  the  core 
K:  intermediate  quantum  number 
TINPO:  the  value  of  TINF  if  CORE-1.5  rr  0.5 

TINF1:  the  value  of  TINF  if  CORE-1.5 

TINF2:  the  value  of  TINF  if  CORE-0.5 

NEFF;  the  effective  principal  o,  antum  number  of  the  outer  electron 

Program  II 

h.  j.  I.  «ULT,  LT,  F.  CORE,  K,  and  NEFF  have  the  same  meaning  as  in 

rogram  I  but  at.  now  arrays.  2,  NEL,  and  C  also  have  the  same 
meaning. 

an  index  running  over  the  states.  In  particular,  it  if  the  index  of 

TUnw.e  StaM  Whe”  <"“lln8  ”Ith  “  tra"Sltl°"-  I"  PROCEDURE 

THOMASIERMI,  however,  it  has  the  value  of  NEFF. 

an  index  running  over  the  final  states. 

in  the  LS  program  only,  means  the  same  thing  as  M  does  in  th.  Jl  and 
IC  programs 


N: 


M: 

K: 
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DNEFF:  the  largest  difference  in  NEFF  for  which  the  program  will  consider 
a  transition  between  two  levels.  Also,  in  the  jl  and  IC  programs, 
Thomas-rermi  wave  functions  are  calculated  for  all  levels  whose  NEFF 
is  less  than  DNEFF&greater  than  NEFF  of  the  ground  state. 

NMAX:  the  highest  energy  level  Index  to  be  considered. 

NNMAX:  the  highest  energy  level  index  to  be  considered  for  the  initial 
state. 

OS:  the  square  of  the  matrix  element  of  r  summed  over  the  Mj  in  the  final 

level  and  averaged  over  the  Mj  in  the  initial  level. 

LAMBDA:  the  wavelength  in  angstroms 

TT:  the  energy  of  a  transition  in  cm  ^ 

Cl:  a  constant  *=0.587234 
YY:  Cl  x  TT 

SUM1:  the  sum  of  the  OS's  for  all  transitions  to  a  given  level 
SUM2:  the  sum  of  the  quantity  OS/TT  for  all  transitions  to  a  given 
level 

SUM3:  the  sum  of  the  quantity  OSxTT  .  r  all  transitions  to  a  given  level 
SI.  S2,  and  S3  have  the  same  meaning  in  tie  LS  program. 

A:  the  transition  rate  in  sec 

II,  12,  and  13  are  indices  used  in  the  Thomas-Fermi  and  intermediate 
coupling  sections  of  tne  jl  and  IC  programs. 

Variables  Associated  with  the  Thomas-Fermi  Method 

QQ:  the  independent  variable  in  the  polynomials  used  to  initialize  the 

ionic  *adius.  See  Stewart  and  Rotenberg, ^Appendix  A. 

QQSUM:  the  sum  of  the  QQ  polynomial 

1/3 

XO:  the  ionic  radius  times  Z  /0.8853 

AA,  BB,  and  CC:  the  coefficients  in  the  initialization  polynomials 

PHIM1 ,  PHIO,  and  PHIP1  are  . ,  and  <J>.  .  respectively.  The  (^-function 

l+l  i  l-i  j/ 2 

is  found  by  solving  the  differential  equation  $''(x)  ■  l<J>(x)/x] 

from  x  inward, 
o 

tf:  the  step  fize  used  in  this  iteration.  It  is  ■  to  X0/4096. 

3/2 

H32:  a  variable  =  to  n 


IISTAN:  the  ste^ire  used  ln  solving  for  (he  actual  futctions_  u 

ywVTC  Vr°  11,6  1“dePendc“t  #*tlabie  i»  the  wave  function  ls 

VMOD=  By2  time,  the  Thomas-Feml  potential  In  Rydbergs.  It  ls  gIve„  as 
a  unction  of  y,  at  points  separated  by  HSTAN/2 . 

WF!  “  2-ilW^1  “tray  giving  R(y2)Ayras  a  function  of  the  Index  of 
energy  leve]  and  y,  at  points  separated  by  HSTAn. 

N0BM:  1t"e  n0r”ll2“£lon  factor  by  which  WF  for  a  given  level  «  Le 
multiplied. 

NTF:  the  index  of  the  highest  energy  level  for  which  a  Thomas-Fermi  wave 
function  is  calculated 

NPOlhTS:  for  each  level,  the  „mber  of  points  at  which  WF  is  given.  It 
is  -  to  the  maximum  y  considered  divided  by  HSTAN. 

NINT:  the  principal  quantum  number  of  the  outer  electron  for  a  given 
state.  (NINT  stands  for  N  integral,  as  opposed  to  N  effective.) 
number  of  nodes  expected  in  the  wave  function  is  NJNT-L-1. 

wri  °  °”ln8  variables  are  declared  insid*  procedure  thomasfermi. 

1,  JFO,  and  WFK1:  values  of  WF  around  a  given  y 
KSP1.  KSO,  and  KSM1:  values  of  Kg  where  the  differential  equation  is 

A  WF(y)  -  K  (y). 

dy  s 

FACTOR?  FACTOAO.  and  FACT0RM1 :  values  of  the  quantity  (1  -  „2r  /12) 
n  and  Ml  stand  for  «  and  -1.  FI  means  that  the  value  is  give! 

or  a  y  one  step  closer  to  the  origin,  while  Ml  corresponds  to  a  y 

one  step  farther  from  the  origin  fTho  «... 

EPSILOv.  8  "•  (The  inte8«tion  proceeds  inward) 

tl'bXLO,.  the  energy  eigenvalue  in  Rydbergs. 

H  (already  declared  outside  the  procedure):  the  step  site  used  in 

solving  for  WF.  It  is  readjusted  as  the  scale  factor  .  is  varied. 

-<2.  a  variable  set  equal  to  H2 
H2T4:  4H2 

H2012:  H?/12 

H2T506:  5H2/6 

Y2T4:  4y2 
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RMAX:  the  largest  value  of  r  considered.  Outside  of  RMAX  the  wave 

function  is  assumed  to  be  too  small  to  make  any  significant  contri¬ 
bution. 

1NVR:  1/RMAX 

SCALF.:  the  scale  factor  a  in  the  scaled  Thomas-Fermi  method 

JUMP:  the  amount  by  which  SCALE  is  varied  between  successive  trials. 

DIRECTION:  the  direction  in  which  SCALE  is  varied.  It  is  +1  for  an 

increase  and  -1  for  a  decrease. 

0LDD1RECT10N:  the  previous  value  of  DIRECTION 

BEFORE:  a  Boolean  variable  which  is  true  as  long  as  SCALE  continues  to 
jump  in  one  direction  by  steps  of  0.05.  When  the  proper  value  of 
SCALE  has  been  overstepped,  BEFORE  becomes  false,  SCALE  converges 
by  a  bisection  procedure. 

LASTTIME:  a  Boolean  variable  which  becomes  true  when  SCALE  has  converged 
sufficiently.  It  activates  the  normalization  and  cut-off  procedures. 

VARYSCALE:  a  Boolean  variable  which  determines  whether  a  scaled  or  un¬ 
sealed  Thomas-Fermi  method  is  to  be  used.  In  the  present  program 
VARYSCALE  is  true  only  for  the  ground  state. 

NEW12:  the  number  of  points  at  which  WF  is  known  after  it  has  been 

interpolated  so  that  the  spacing  between  points  is  HSTAN 

1.P0LY:  a  function  of  L  which  is  used  a  large  number  of  times.  It  is  = 
to  16L(L+l)+3. 

NODES:  the  number  of  times  that  the  wave  function  crosses  the  axis 

each  time  the  differential  equation  is  solved. 

RATIO:  the  ratio  of  the  values  of  the  wave  function  between  the  two 

points  closest  to  the  origin. 

DRATIO:  tne  change  in  RATIO  between  two  successive  iterations. 

INSIDE:  a  variable  which — on  the  last  iteration — goes  to  true  as  soon 
as  y  moves  from  the  exponential  region  into  the  oscillatory  region. 

N0RM0D:  a  non-array  quantity  which  *•  summed  to  form  the  array  NORM, 
the  normalization  factor. 
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Variables  Associated  with  Intermediate  Coup] tng 

(Not  a]l  of  the  variables  declared  are  described  here,  but  most  of  the 

important  ones  are  mentioned.) 

SPORBMATRIX  (4D  or  3D):  the  coefficients  of  the  first  spin-orbit  para¬ 
meter  in  the  theoretical  Hamiltonian  matrix.  In  the  4D  case,  the 
first  index  is  the  1  of  the  configuration,  the  second  index  is  the 
J,  and  the  last  two  are  the  row  and  the  column  within  a  J-block. 

The  3D  (i.e.  3-dimensional)  array  is  missing  the  first  index. 

SPORBMATRIXP  (4D  or  3D):  the  coefficients  of  the  second  spin-orbit 
parameter  in  the  Hamiltonian. 

ELECMATRIX  (4D  or  3D):  the  theoretical  Hamiltonian  for  the  electro¬ 
static  parameters.  In  the  4D  array,  the  first  three  indices  are 
the  same  as  in  SP0RBMATR1X4D  while  the  last  index  identifies  the 
electrostatic  parameter  under  cons ' deration.  No  -olumn  *-Hex  is 
necessary  since  the  electrostatic  Interaction  is  diagonal  m  the 
LS  basis. 

IELEC  and  NELEC:  the  index  and  total  number  of  the  electrostatic 
parameters. 

IPARAM  and  NPARAM:  the  index  and  total  number  of  parameters. 

ELEC.PARAM:  an  array  containing  the  values  of  the  electrostatic  para¬ 
meters  as  they  are  adjusted  by  PROCEDURE  DIAG. 

SPORBPARAM  and  SPuKBPARAMP:  the  values  of  the  two  spin-orbit  parameters 

I GREATER :  the  greater  of  the  two  £'s  of  two  states  under  consideration. 

STLS3D:  the  matrix  elements  of  cos6  between  pairs  of  LS  states.  The 

trrac  index  is  LGREATER  and  the  other  two  are  the  rows  and  columns. 

J2,  J3,  S2,  S3,  L2,  and  L3:  the  variables  used  to  represent  quantum 

numbers  in  generating  the  Hamiltonian  matrices  and  the  LS  strengths 

HOLD:  the  variable  that  holds  on  to  the  value  of  L  as  the  energy  levels 
are  read  in.  When  the  end  of  a  configuration  is  reached,  the  next 
L  is  no  longer  ■  to  HOLD  and  PROCEDURE  DIAG  is  called. 

T2D  and  N2D:  arrays  holding  the  values  of  T  and  N  for  the  levels  of  a 
given  configuration.  The  first  index  is  the  J  of  the  level  and  the 
second  puts  the  states  of  a  given  J-block  in  order  of  increasing 
energy . 
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IJ,  IJMIN,  and  UMAX:  the  index  and  limits  for  the  J-blocks  of  a  given 
configuration. 

INDEX  and  BLOCK:  arrays  over  N  giving  the  J-block  and  the  position  in 
that  block  of  every  level  fed  into  the  program.  The  number  BLOCK 
is  distinct  from  one  configuration  to  the  next. 

LOWERBLOCK:  the  lowest  BLOCK  in  a  given  configuration. 

SIZE2  and  SIZE3:  the  sizes  of  the  J-blocks  as  predicted  by  the  sections 
which  generate  the  theoretical  Hamiltonian  matrices  and  the  LS 
strengths . 

SIZE:  the  size  of  the  J-blocks  in  the  observed  levels  which  are  read  by 

the  program. 

I2L0WF.R  and  I3L0WER:  integer  arrays  which  indicate  what  position  in  the 
STLS3D  matrix  to  look  for  a  given  pair  of  L's  and  J's. 

TRACE  and  TRACES:  the  traces  of  the  J-blocks.  Since  the  traces  of  the 

undiagonalized  theoretical  Hamiltonian  and  the  traces  of  the  observed 
energy  matrix  should  be  the  same,  enough  information  is  present  to 
initialize  the  values  of  the  parameters  for  a  given  configuration. 

H:  a  matrix  for  the  Hamiltonian  formed  from  the  theoretical  expressions 

for  the  Hamiltonian  and  the  values  which  the  parameters  happen  to 
have  t.t  the  moment. 

TCALC:  :n  array  containing  the  energies  calculated  by  diagonalizing  the 
H-matrix. 

V3D,  V2D,  and  VTEMP:  various  forms  of  the  matrix  which  diagonalizes  the 
Hami ltonian. 

DERIV:  a  matrix  of  the  derivatives  of  the  calculated  energy  eigenvalues 
with  respect  to  the  parameters. 

ALPHS:  the  coefficients  of  the  set  of  simultaneous  equations  which  are 

solved  for  the  corrections  to  the  parameters.  They  are  calculated 
from  the  DERIV' s  and  the  differences  between  the  observed  and 
calculated  T's. 

IT:  the  number  of  times  that  the  parameters  for  a  configuration  have 

been  readjusted  for  a  better  least-squares  fit  to  the  observed 
energies . 
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3JMSQUARES:  f  -  sum  of  the  squares  of  the  differences  of  the  observed 
and  calcuxated  energies. 

HCLDSUM:  the  value  of  SUMSQUARES  from  the  previous  iteration.  If 

S UNSQUARES  has  changed  by  less  than  0.1%,  the  iteration  is  terminated. 


REAL  PROCEDURE  XI  is  the  same  as  -(J,  J’) 


(J  +  1)(2J  +  3)  if  J'  =  J  +  1 

=  J(J  +1)  J*  «  J 

J(2J  -  1)  J'  =  J  -  1 

\ 

REAL  PROCEDURE  W2MOD  is  a  modified  Racah  coefficient  squared. 

W2M0D(L,  J,  L\  J',  S)  -  (2J '  +  1)  L  (21.  +  1)(2L’  +  1)  W2  (L,  J,  L',  J';S) 
Moreover, 


W2M0DC1.  J.  L\  J1.  S)  I <SLJ  P  SL ' J  1  >  1  2 


XI(J,  J’) 


<L  P 


L  *  >  I  2 


where  P  is  class  T  with  respect  to  J  and  [S,  P]  *  0 


REAL  PROCEDURE  I  calculates  the  radial  integral  J  R  (r)R  (r)  rdr  by 

numerical  integration  it  both  N  and  M  are  le3S;  than  the  Jndex  of 
the  highest  state  for  which  a  Thomas-Fermi  wave  function  has  been 
calculated.  For  all  other  pairs  the  Bates-Damgaard  method  of 
integration  is  called  through  REAL  PROCEDURE  BD. 

PROCEDURE  THOMASFERMI  takes  NEFF,  NINT,  and  L  for  the  Kth  level  and 
produces  a  wave  function. 

PHI/GENERATOR  produces  the  Thomas-Fermi  potential  by  numerically  solving 
a  differential  equation  derived  from  the  theory.  All  quantities 
calculated  in  this  section  perta.in  to  all  levels  for  which 
PROCEDURE  THOMASFERMI  is  applied. 

Section  in  the  program  which  sets  NINTtN]  for  all  the  energy 

levels  for  which  Thomas-Fermi  wave  functions  are  calculated 
must  be  changed  for  every  different  atomic  or  ionic  species. 
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The  IC  Program 

MaTRIX/GENERATOR  -  the  part  of  program  which  generates  the  theoretical 
expressions  for  electrostatic  and  spin-orbit  perturbations  to  the 
Hamiltonian.  Th?  LS  basis  Is  used  and  matrices  similar  to  those 
on  pages  268-269  of  Condon  and  Shcrtley  i11  are  produced  for  each 
parameter.  Section  1  3  page  299  gives  the  actual  coefficients 

generated,  since  the  program  is  at  present  designed  to  do  the  rare 
gases. 

LS/ STRENGTHS  -  produces  a  matrix  for  the  angular  part  of  the  dipole 
operator.  All  possible  pairs  of  LS  states  are  considered.  When 
PROCEDURE  DIAG  has  expressed  each  energy  level  as  a  lines-  combina¬ 
tion  of  LS  states,  it  is  then  easy  to  calculate  the  angular  matrix 
element  between  any  two  of  the  actual  states. 

A  section  .n  the  program  takes  the  energy  levels  in  each  configuration 

as  a  group,  places  them  in  blocks  of  different  J,  and  then  re-orders 
each  J-block  in  order  of  increasing  energy,  which  are  then  fed  into 
PkOCEDURE  i  IAG,  ore  configuration  at  a  time. 

PROCEDURE  DIAG  takes  each  configuration  fed  into  it  and  compares  the 
actual  obserad  energies  with  the  theoretical  expressions  for  the 
Hamiltonian  which  have  been  generated  by  HATRIX/GENERATOR .  DIAG 
then  adjusts  the  theoretical  parameters  so  that  the  diagonalization 
of  the  theoretical  Hamiltonian  gives  a  least-squares  fit  to  the 
observed  energy  igenvalues.  Once  this  is  done,  the  matrix  which 
diagonalizes  the  Hamiltonian  is  merely  an  expansion  of  tha  observed 
levels  in  of  LS  basis  states.  This  matrix  saved  for  the  final 

computation  of  the  oscillator  strengths. 
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INPUTS  TO  PROGRAMS  1  AND  2  OF  GIANT/ FIASCO  (LS  PP  GRAM) 


under 


under 


"FIASCO"  (input  to  Program  i): 

a  card  with  Z,  the  atonic  number  and  N  ,,  the  number  of  electrons 

el 

a  card  with  the  value  of  T,  , 

inf 

for  each  level  of  Moore's  tables  a  card  giving; 

L  (orbital  angular  momentum  of  the  outer  electron) 

MULT  (multiplicity  =  2S  +  1) 

LT  (total  orbital  angular  momentum) 

P  (parity  +  +1  even,  *  -1  odd) 

J  (total  angular  momentum) 

T  (energy  of  the  level  in  cm  ^) 

"READER2"  (input  to  Program  2); 

a  card  with  NMAX  (the  maximum  N  to  be  read)  and  NNMAX  (max  on  the 
index  of  the  initial  state) 

a  carJ  with  DNEFF  (the  limitation  on  interacting  levels) 


INPUTS  TO  PROGRAMS  1  AND  2  OF  GIANT/FIASCO  (jl  and  IC  PROGRAM) 


under 


under 


"FIASCO"  (input  to  Program  D • 

a  card  with  Z,  :he  atomic  num.  . c  and  N  ,,  the  number  of  electrons 

el 

a  card  with  the  value  of  T,  ,  for  CORE  *  0,  the  value  of 

inf 

T,  ,  for  CORE  *  1.5,  the  value  of  T,  -  for  CORE  *  0.5 
inf  inf 

for  each  level  of  Moore's  tables  a  card  giving: 

CORE  (total  angular  momentum  of  the  core) 

L  (orbital  angular  momentum  of  the  outer  electron) 

K  (intermediate  quantum  number) 

J  (total  angular  momentum) 

T  (energy  of  the  level  in  cm  ^) 

"READERS"  (input  to  Program  2): 
a  card  with  NMAX  (the  maximum  N  to  be  read) 
a  card  with  DNEFF  (the  limitation  on  interacting  levels) 
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15  cards  giving  the  co.fflcl.nt.  rf  the  three  lnitiallratlon 

polynomials  for  the  Thomas-Feml  function.  Each  card  has 
four  numbers: 


the  index  of  the  coefficient, 
A (Index) , 

B (index) , 


index) , 

(See  Stewart  and  Rotenberg5,  Appendix  A) 

and*  &-&£-!£-  program  onlg,  for  each  possible  L  (L  *  0,  1,  2  3) 

a  card  with: 

L  -  1), 

L), 

L  +  1), 


LMj.  (coeff  of  F„  if  L 

2  total 

LSAME  (coeff  of  F  <f  l 

2  *  total 


LP1  (coeff  of  F0  if  L 

2  total 


LM1SC  (coeff  of  u  if  I  _  T  , 

u0  if  Ltotal  L  1  and  S  =  0), 

LP1S0  (coeff  of  G2  if  Ltot  ^  “  L  +  1  and  S 


-  LUl  L 

VS>ee  Condon  and  Shortley,  pp  2SJ-299). 


Cx 

-  t 


rnuuKAFld 


is  Program:  For  elements  In  which  the  L  of  the  parent  tern  la  non-rero. 

t  e  expression  for  OS  must  be  changed  fro.  W2M0D  (1,  J,  L-,  j.(  s)„2  u 


"2M0D  (L,  J,  L'f  J  ,  S) 


W2M0D  (1,  L,  1',  L'  L  ) 

- - - parent  2 

XI(L,  L’)  0 


In  more  complex  element.,  coefficients  of  fractional  parentage  must 
also  be  added. 

J1  and  IC  Programs:  To  do  other  rare  gases,  the  only  change  necessary 
IS  the  setting  of  HINT  (the  principal  quantum  number)  for  all  of 
the  levels  to  be  done  by  Thomas-Fermi  method. 

To  do  other  elements  the  cards  which  define  the  ground  state  as  a 

Sc  state  and  multiply  the  strengths  of  ground  state  transitions  by 
six,  must  be  altered. 
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IC  Program:  NMAX  must  be  set  low  enough  that  no  configuration  with  fewer 
than  the  predicted  number  of  levels  gets  read  by  Program  2.  NM.'X 
should  be  one  greater  than  the  last  le-  j1  in  the  highest  configura¬ 
tion  to  be  done.  Although  this  grea  ly  limits  the  IC  Program,  it 
would  be  possible  to  do  many  more  transitions  if  PROCEDURE  DIAG 
expanded  the  actual  states  in  terms  of  jl  basis  states  instead  of 
LS  states.  Higher  states  in  incomplete  configurations  could  then 
be  approximated  as  pure  jl  states.  Therefore,  a  matrix  which 
transforms  the  LS  basis  to  the  j basis  would  be  a  useful  addition. 

For  elements  with  more  complex  configurations,  the  viriable  which 
identify  the  states  and  the  J-blocks  must  be  changed,  alon  ,ith 
tb  ■  generator  of  the  theoretical  Hamiltonian  and  the  LS  st  engtht.. 
This  involves  rewriting  most  of  the  procedures. 
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OUTPUT  FROM  PROGRAM  2  OF  GIANT/FIASCO 


Following  the  listing  of  the  energy  levels  from  Program  1, 

Program  2  will  print  out  varying  quantities  of  fairly  self-explanatory 
information,  depending  upon  the  coupling  scheme — that  is,  LS,  jl,  or  IC. 

The  jl  and  IC  programs  will  first  list  all  the  Thoma ''-Fermi 
wave  functions  s s  a  function  of  the  number  of  steps  from  the  origin  in 
y.  (The  function  printed  is  actually  q(y);  see  the  section  describing 
GIANT/ FIASCO. )  Note  that  the  first  ten  steps  are  given,  and  then  jumps 
of  five  steps  are  taken  to  save  space.  Immediately  below  each  wave 
function  is  printed  the  expected  number  of  nodes,  which  is  simply  n^n  -1-1, 
and  the  actual  number  of  nodes  in  the  wave  function  shown.  The  two 
numbers  are  always  the  same  for  levels  done  by  the  scaled  TF  method,  but 
in  the  unsealed  cases  (These  can  be  identified  by  the  fact  that  SCALE  = 
1.000.)  there  is  no  guarantee  that  this  is  so.  When  the  two  are  the 
same,  it  is  an  indication  that  the  Thomas-Fermi  method  is  valid. 

The  IC  program  will  now  print  out  the  transformation  matrices 
from  LS  states  to  the  observed  states.  Double  lines  separate  different 
configurations,  wh^’e  single  lines  separate  the  iterations  approaching 
a  leastsquares  fit  for  a  given  configuration.  For  each  block  of  a  dif¬ 
ferent  J,  the  observed  energies  and  the  energies  calculated  by  diagonal¬ 
izing  the  theoretical  Hamiltonian  appear  next  to  the  corresponding  rows 
of  the  transformation  matrix.  (Note  that  some  of  the  matrices  are 
rectangular  although  all  of  them  should  be  square.  Extra  columns  on 
the  right  should  be  ignorea.)  Directly  below,  and  within  the  same  J 
block,  are  the  derivatives  of  the  calculated  energy  eigenvalues  with 
respect  to  the  various  adjustable  parameters.  The  rows  are  in  order  of 
increasing  energy  eigenvalues,  and  the  columns  are  in  the  same  order  as 
the  list  of  parameters,  which  appears  once  for  each  iteration.  (Deriva¬ 
tives  for  parameters  which  are  not  involved — that  is,  parameters  which 
are  zero  after  the  first  iteration — should  be  ignored.) 

Finally,  the  LS,  Jl,  and  IC  programs  all  print  out  the  transi- 

2 

tion  probabilities  in  the  form  3 1 < n J  z | m> |  ,  which  had  already  been  summed 
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over  final  states  for  m  and  averaged  over  initial  states  for  n.  At  the 

far  right  is  Einstein’s  spontaneous  emission  coefficient  in  the  fom 
—8 

x  10  .  However,  at  present  the  jl  and  IC  programs  have  a  card  near 

che  end  of  Program  2  which  reads: 

A-*-Ax(2xJ (N)+l/  (2xJ (M)+l) ; 

_8 

This  reverses  the  indixes  n  and  m  so  that  A  x  10  is  printed  out 
-8  mn 

instead  of  A  x  10  as  the  comment  cards  in  Program  2  claim, 
nm  6 

AUXILIARY  PROGRAM  FIASCO /RESULTS 

This  pro$  ram  merely  reads  the  stored  output  of  Program  2  from 
the  disc  READF2  and  prints  out  the  transition  probabilities  in  several 
forms.  GF  is  gf,  FNM  is  f  ,  ANM  is  A  ,  and  AMN  is  A  .  Input  is  just 
the  levels  of  Moore's  tables  with  the  same  cards  that  go  into  Program  1. 
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APPENDIX  C 


CALCULATION  OF  THE  OPTICAL  RADIATION  EMITTED  BY  A 
CYLINDRICAL  AFC  OF  KNOWN  TEMPERATURE  PROFILE 

by.  R.  Liebermann 


INTRODUCTION 

In  order  to  compare  measurable  physical  properties  of  an 
emitting  plasma  with  theory,  a  computer  program  has  ben  written  to 
calculate  at  the  surface  of  a  cylindrical  arc  discharge: 

(1)  the  spectral  radiance  (normal  to  the  surface)  due  to 
radiant  contribution  from  along  the  line  of  sight  of 
the  diameter  which  has  a  spatial  temperature  gradient, 
and 

(2)  the  spectral  radiant  emittance  due  to  contributions  from 
all  radiating  volume  elements  of  a  homogeneous  temperature 
arc  column. 

The  above  calculations  involve  self-absorption  processes  within  the  arc 
column  and  therefore  are.  solutions  of  a  radiant  heat  transfer  equation. 


A.  Radiant  Heat  Transfer  Equation 

With  the  assumption  that  the  plasma  is  in  local  thermal 
equilibrium  (LTE)  the  spectral  radiance  is  related  to  the  properties  of 
the  medium  through  the  radiant  heat  transfer  equation 


-  d  I 
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d  r 
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(1) 


which  describes  the  spatial  gradient  of  intensity  along  the  line  of 
sight  at  some  point  x  where  1^  is  the  spectral  radiant  intensity  at 
frequency  v  in  the  proper  direction,  is  the  corresponding  black 
body  intensity  and  is  the  optical  depth  as  measured  from  x  =  -n0  to  x. 
We  define  the  r\  tical  depth  as 
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»rV”)  T  18  thC  SP":tral  abSOr!,tiVit>’  (including  stimulated 

emission)  due  to  free-free,  bound-free,  end  bound-bound  radiative 

processes.  The  solution  of  1  ib  given  as 


where 
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(3) 


(4) 


with  L  corresponding  to  that  location  along  t.e  line  of  sight  (at  the 
surface  of  the  arc  column)  having  spectral  radiance  I  (t  )  .  In  3  it  was 
assumed  for  x  <  0  that  no  contribution  to  I  (x  )  is  made  (i.e.,  no 
external  radiation  field  is  present). 

1 •  Spectral  Radiance 

If  our  interest  is  radiance  at  the  arc  column's  surface  due 
to  contributions  from  along  the  line  of  sight  of  the  diameter,  equation 

by  an  assumption  of  radial  symmetry  can  be  placed  in  the  following 
form.  6 


VV  •  Vr-R>  - 


By  (Z) 
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exp  {  Z  }  +  exp  {-Z} 


exp  i -Z )  dZ 


(5) 
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where  Z  is  an  optical  depth  measured  from  the  center  of  the  arc  column 
to  some  radial  point  r  along  the  radius  which  bar,  a  maximum  dimension 
R.  The  correspondence  is  Z  »  o  at  r  *  o  and  Z  =  2  at  r  =  R  thus 
fR  rr 


Z  = 


(T(r) )  dr  ;  Z  = 


*1  (T(r) )  dr 


(ft) 


For  programming  purposes  Equation  5  is  the  better  form  of  3,  In  order 
to  solve  5  it  is  first  necessary  to  have  prior  knowledge  of  the  tempera¬ 
ture  profile  T(r)  and  for  a  given  gas  vapor  pressure  a  method  by  which 

K  '  can  be  determined  such  that  solutions  of  6  are  obtainable.  Part  D 
v 

of  this  section  describes  the  general  theory  involved  for  determining 

< y(T). 

Utilizing  the  theory  of  Part  D  we  are  then  prepared  to 
numerically  integrate  Equations  6  and  thus  5.  The  method  of  approach  was 
as  follows.  For  a  specified  temperature  profile  and  frequency  the 

integrand  of  6  was  determined  at  division  points  in  r  and  interpolated 

2  2 
in  r  between  division  points.  Using  a  linear  expansion  of  < ^  in  r 

a  closed  form  integration  was  performed  between  division  points  and  a 

summation  yielded  Z(r).  In  similar  fashion  the  black  body  term  in  the 

integrand  of  5  was  calculated  at  divisioi  points  in  r  and  linearly 
2 

expanded  in  Z  .  A  closed  form  integration  of  5  was  then  possible 
between  division  points  of  Z(r)  and  was  summed  to  yield  T^ft^).  Although 
the  latter  integration  can  be  performed  for  all  values  of  the  parameters 
serious  errors  (due  to  cancellation  effects)  may  develop  when  AZ  is 
small  between  end  points  of  the  intervals.  Hence  for  AZ  <.01  the 
solution  was  expanded  in  terms  of  AZ  whereas  for  AZ>.01  the  closed 
form  was  used.  Solution  of  5  was  then  repeated  for  other  selected  v 
values  of  the  specified  band  pass  to  yield  a  representative  spectrum. 

2 .  Spectral  Radiant  Emittance 

If  our  interest  is  spectral  emittance  at  the  arc  column's 
surface  due  to  contributions  from  all  radiating  volume  elements  with  the 
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discharge  assumed  homogeneous  in  temperature  then  Equation  3 
which  yields  the  spectral  radiance  along  a  particular  line  of  sight 
reduces  to 


IV(T,S) 


■„<T> 


1  -  exp  { 


where  now  from  Equation  4 


(7) 


t  -  K'  (T)  S 
v  v 


(8) 


S  being  the  distance  through  which  radiation  travels  along  a  particular 
line  of  sight  and  is  a  function  of  arc  geometry. 

The  spectral  radiant  emittance  at  the  surface  of  the  arc 
column  can  be  calculated  using  7  by  applying  the  cosine  law  and 
integrating  over  solid  angle.  Using  the  geometry  shown  in  Reference  1 
the  spectral  radiant  emittance  is  given  as 


Fv(T) 


VT> 


■ir/2 

1  r 

v  -  4 

exp  /  — r — ^ — r - r~  >  udud$ 

J 

3 

!  u  f  (1-u  )  Sin  J 

(9) 


where  D  is  the  column  diameter 

The  radiant  emittance  is  then  determined  by 


F(T) 


>oo 


Fy(T)  dv 


(10) 


Equations  (9)  and  (10)  were  solved  using  straightforward  numerical 
techniques. 
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B.  Input  to  Arc  Intensity  Program 


The  Input  required  for  the  above  calculations  is  fully 
described  within  the  program  itself  which  is  shown  in  Part  E  of  this 
Section.  In  brief  the  basic  required  input  is  as  follows: 

(1)  Spectroscopic  data  involving  bound-bound  transitions 
such  as 

(a)  oscillator  strengths 

(b)  electron  impact  line  shifts  and  widths 

(c)  excitation  energies 

(d)  statistical  weights 

(2)  thermodynamic  properties  uuch  as  tables  of 

(a)  electron  number  densities 

(b)  neutral  particle  number  densities 

(c)  internal  partition  functions 

(d)  ionization  lowering  values 

(3)  temperature  profile 

(4)  specified  wavelength  band  pass. 

C.  Output  of  Arc  Intensity  Program 

The  output  of  the  program  provides  plots  as  well  as  printouts 
of  the  spectral  radiance  and  spectral  radiant  emittance.  In  addition 
plots  of  spectral  absorptivities  are  obtainable.  Samples  of  these 
plots  are  shown  in  Figures  3,  4 ,  and  5.  For  particulars  on  how 
to  obtain  these  plots  see  program  Part  E. 

D.  Spectral  Absorptivities 

The  net  absorption  coefficient  which  is  the  difference 
between  true  absorption  and  induced  emission  coefficients  is  in  a 
spectral  line  for  a  given  frequency  v  and  temperature  T  (assuming  LTE) 

'  *v  N„  f1-""  pr}]L(v)  <u) 
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where  it,  r^,  c,  h  and  k  are  the  usual  constants,  fmn  is  the  absorption 
oscillator  strength  for  transitions  from  state  n  to  m  (n  being  the 
lower  energy  level  state),  is  the  population  density  of  the  initial 
state,  and  L(v)  is  the  line  profile  normalized  such  that 


L(v)dv 


(12) 


where  L(v)  is  a  function  of  the  lines,  shifts,  widths  etc.  all  a  function 
of  temperature  and  pressure. 

For  conditions  of  LTE  the  population  density  can  be 
expressed  in  terms  of  $o  (the  total  number  of  atoms  or  ions  in  question) 
by  use  of  the  Boltzmann  relation 


g  N 
n  o 


exp 


(13) 


where  g  is  the  statistical  weight  of  the  initial  level,  E  is  the 
n  n 

excitation  energy  of  the  initial  level, Uq  iR  the  atomic  (or  ionic) 
internal  partition  function  and  T  is  the  absolute  electron  temperature. 

The  line  profile  L(v)  is  course  dependent  on  the  mechanisms  of 
line  broadening.  The  stark  effect  due  to  electron  impact  has  the 
Lorentz  dispersion  profile  given  as 


L(v)  “  ± 
HU) 


1 _ ' 

fv-(v  +d)  ^  1 

I  n 


(14) 


where  u>  ia  the  half  width  of  the  line  at  half  intensity  and  d  is  th<_ 
electron  impact  shift,  is  the  observed  frequency  of  the  line  in  an 
uuperturbed  condition. 

When  resonance  broadening  mechanisms  are  of  importance,  one 
must  add  to  the  electron  impact  half  width  u>  the  resonance  half  width 
ui  which  from  first  order  approximation  can  be  expressed  for  resonance 

lines  as 
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\  1/2 


where  g  and  g^  are  respectively  the  lower  and  upper  states'  statistical 
weight  factors  of  the  resonance  line,  N  is  the  ground  state  population 
density  and  f  is  the  absorption  oscillator  strength  of  the  transition 
(actually  f  as  defined  here  is  the  average  value  for  the  resonance 
multiplet  transition  obtained  by  a  sum  of  the  multiplet's  individual 
f's  over  final  states  which  are  then  averaged  over  initial  states).  In 

equation  15,  N  are  neutral  perturbers  of  the  same  kind  (radiating  atoms). 

8 

For  the  case  when  perturbers  of  unlike,  kind  are  of  such  abundance 
(e.g.  a  strong  concentration  of  another  gas  specie  is  present)  then 
Van  der  Waals  broadening  may  be  of  some  importance  and  should  be 
included  in  the  dispersion  profile. 

Besides  line  broadening  due  to  charge  perturbers  (Stark 
effect)  or  neutral  atoms  (Resonance  and  Van  der  Waals  effects)  there 
exists  other  broadening  mechanisms  which  must  also  be  considered, 
namely  Doppler  broadening,  which  results  in  a  Gaussian  line  shape 
profile.  This  type  of  profile  is  dominant  when  collisions  are  negligible 
and  the  velocities  of  the  radiating  systems  have  a  thermal  distribution. 
The  Gaussian  line  shape  is  given  as 


ld<^>  -  prj 


—  exp 

to 
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-?n2  t A 
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with  Av  =  v— v  and  where  ui  is  the  Doppler  half  (half)  width 
o  D 
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2kTln2 

i)  m  ' 

D  Mc2 

m 


M  being  the  mass  of  the  radiating  system.  When  collision  effects  cannot 

be  neglected  (e.g.  when  w,  u>  etc.  are  of  the  same  order  of  magnitude  as 

K 

Wn)  then  joint  effects  of  Doppler  and  the  more  important  pressure 
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broadening  mechanisms  should  be  considered.  The  observed  line  shape 
is  then  the  so  called  Voigt  profile  provided  the  broadening  processes 
are  statistically  independent.  From  17  it  is  obvious  that  lines 

occurring  in  the  Uy  regions  will  tend  to  be  more  of  a  Voigt  profile 
type. 


The  net  absorption  coefficient  k/  as  expressed  in  .6  is  not 
only  due  to  bound-bound  transitions  but  also  includes  bound-free  and 
free-free  processes  as  well.  Therefore,  one  must  add  to  a  sum  of 
Equation  11  these  other  processes  to  obtain  representative  values  of 
5  in  the  spectral  regions  wheie  bound-bound  transitions  are  not 
dominant.  For  conditions  of  LTE  the  bound  free  absorptivity  K^f(vtT) 
can  be  expressed  in  terms  of  a  sum  over  excited  states  designated  by 
quantum  numbers  n  and  l  having  corresponding  photoionization  cross- 

sections  anl  and  population  density  N  whose  edge  frequencies  v 

i  ,  ,  “»*  n£ 

are  less  than  or  equal  to  v 


(18) 


The  net  bound-free  absorptivity  is  then  expressed  as 


(19) 


Solution  of  18  is  a  formidable  task  since  for  each  level  requires 
numerical  solutions  of  the  Schrodinger  equation.*  A  simplified  approach 
is  to  separate  18  into  two  distinctive  regions  i.e,,  into  levels  which 
are  hyurogen  like  and  into  levels  which  are  non  hydrogen  like  thus 


bf 


(v.T) 


(20) 


* 

See  foi  example  Ref.  1, 
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where  in  general  for  the  hydrogen  like  terms  n>>l.  The  sum  over  the 
hydrogen i<:  levels  can  be  simplified  further  as  follows.  Starting  with 
a  sum  rule  analogous  to  that  for  bound-bound  transitions 


H, 

o  .dv 
n£ 


irr  cf  . 
o  n£ 


H  -.3 

and  using  the  proportionality  o  '  v  the  hydrogenic  photo ionization 

H  n*’ 

cross  section  0  ,  is 

nil 

H  Vn92 

oH0  -  2nr  c  f  .  (22) 

nil  o  3  nt 

v 

where  f  corresponds  to  the  excited  n,£  state's  continuum  oscillator 
strength.  Using  the  Rydberg  formula  for  hydrogenic  energy  levels,  the 
binding  energy  E„  of  the  excited  level  is  (E  .  I  E  ) 


where  Z  refers  to  the  ionization  stage,  n  is  the  principal  quantum  number, 
E^  is  the  ionization  potential,  E  is  the  excitation  energy  of  the 
excited  level  and  E  is  the  Rydberg  constant.  From  23  we  see  that  the 
rate  at  which  the  binding  energy  changes  with  respect  to  unit  change  in 
principal  quantum  number  for  n>>l  is 


aeb  - 


-2ehz 


We  new  define  an  average  hydrogenic  cross  section  for  all  £  states 
with  same  quantum  number  n  as 


(2£+l) 


(25) 
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(2£+l)  o 
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where  £  is  the  orbital-angular  momentum  quantum  number .  Thus  using 
the  hydrogenic  relations 


and 


I 

9 

I 


(2£-rl)f  :  .5n 

nt 


(2£+l) 


n=constant 


(26) 


Equation  25  takes  the  form  '  f 


o 

•i 


o 

4H*- 


(27) 


where  a  and  aQ  are  the  usual  constants.  Equation  25  is  the  familiar 

hydrogenic  cross  section  less  the  Gaunt  factor  G,,.  Using  Equation  13, 

br 

and  not.'ng  for  atomic  gases,  who  e  higher  levels  are  hydrogen  like  but 
degenerate  in  J  (total  angular  momentum  quantum  number)  that 


8 


n 


(2J+1)n 


(28) 


where  is  the  statistical  weight  factor  of  the  ground  state  of  the 
parent  term;  and  using  the  relationship  of  Equation  24  the  hydrogenic 
sum  term  of  20  becomes 


k  |V  n  y  j-y 

UQ  Eu  Ihv)  °Z^eXp\  kT  j 


&Er 


•  G 


bf 


(29) 
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For  hydrogen  the  excitation  energy  is  approximately  equivalent  to 
E^.  However,  for  other  gases  whose  higher  levels  are  hydrogen  like 
but  degenerate  in  J  then  should  be  averaged  over  these  J,£  states 
which  have  E  excitation  energies  such  as 

E  -  )  (2J+1)  E  /2n2g  (30) 

n  /  n  J ».  i 

J 


Thus  with  E  as  defined  in  30  obeying  Equation  23  wherein  E 

n  a  l 

replaced  by  E  then  29  takes  the  form  of 


L—> 

n  l 


H  H 

o  „  N  „ 
n  H  at 


/  2  2 
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U  E.. 

o  1 


i  E  13 

-h!l  Z\Cbf 


f_EJ! 

exp 


is 


\ 

exp 


(31) 


with  the  assumption  that  =  constant  for n  >>1. 

The  sum  on  the  right  of  31  may  for  practical  purposes  be  taken  as  an 

integral  with  lower  limit  E  *  hv  and  upper  limit  E  -o.  Hence  for 

on  d 

representing  the  maximum  edge  frequency  of  the  merged  levels  where 
n>>l  with  v  <  v  then 


or  for  the  case  when  v  is  less  then  the  maximum  edge  frequency  of  the 
merged  levels  (v<v  )  then 


C-ll 
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'tU  nZ 


1.82  x  10 


24 


111 
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1  kTJ  rp 


jbf 


(33) 


Eq..  tlcns  32  and  33  are  equivalent  to  those  presented  by  Biberman, 
Norman  ar.d  Ulyano^\ess  their  correction  factor  which  for  cur  high  n 
levels  is  approximately  one  and  equivalent  to  the  Gaunt  factor.  It  is 
of  interest  to  note  that  in  the  form  of  31  if  the  right  hand  sum  is 
taken  as  an  integral  with  lower  limit  set  to  E  -o  and  upper  limit 

D 

E  ■  -  •  representing  unbound  electrons  we  obtain  at  once  the  familiar 

D 

free-free  absorptivity  expression 


'ff 


1.82  j:  10 
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(34) 


where  the  term  has  been  introduced  and 
Thus  equation  20  takes  the  form  (including 


Kbf+tf 


,  a  N  „  +  1.82  s  10 
/  .  r.  n£ 


24 


is  the  free-free  Gaunt  factor, 
the  f-f  absorption) 

2  r 


gi T  V 


U  V 
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exp 


;-e~ 


kT 
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(  exp 


(35) 


whereby  previous  definitions  v'  =  v  if  v>  v  and  v'  =  v  if  v<v  . 

n  ”  n  n 

In  equation  35  the  non-hydroganic  sum  term  then  only  includes  such 

levels  having  edge  frequencies  v  which  meet  the  requirement  that 

v  <  v_  i  v.  In  the  above  derivations  the  assumption  that  Grc  is 
n  E  ff 

constant  independent  on  temperature  and  frequency  is  valid  provided 
T  <  2  ev  (see  Ref.  3  )  however  when  T  is  of  greater  value  than  G^ 
increases  and  varies  somewhat  with  frequency.  Likewise  when  n  is 
of  some  lower  value  then  G^  is  no  longer  a  constant  but  varies  with 
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quantum  number  and  frequency.  The  accuracy  of  35  is  thus  better 

enhanced  when  the  sum  term  includes  the  majority  of  the  excited  levels. 

Such  calculations  involving  the  sum  term  have  been  made  by  Ref.  2 

and  ^  who  have  utilized  the  techniques  proposed  by  Ref.  1  for 

determining  photoionizat ion  cross  sections  for  ionized  gases  other 

than  hydrogens.  However  their  results  differ  somewhat  in  that  they 

have  chosen  different  normalization  factors.  In  essence  (using  their 

results)  Equation  35  remains  the  same  with  the  exception  of  replacing 
2 

the  factor  Z  G,  .  with  a  t  (Zeta)  factor,  and  in  addition  when  using  the 
bt 

t;  factor  of  Ref.  4  .  v'  is  to  be  replaced  by  v  for  all  values.  Their 
sum  term  then  only  includes  a  limited  number  of  lower  excited  levels 
(including  the  ground  state). 

Since  several  b-b  transitions  usually  occur  for  a  given  lower 
excited  state  near  the  frequency  f  its  continuum  photoionization 
edge  (the  upper  excited  state  being  near  the  ionization  limit) 
there  results  a  pseudo  continuum  (i.e.  the  lines  in  question  merge  with 
the  continuum  causing  in  appearance  a  false  photoionization  edge  at 
a  slightly  lower  frequency).  When  these  b-b  transitions  are  not  included  in 
the  sum  of  Equation  11  near  such  corresponding  photoionization  edges 
then  serious  errors  may  develop  in  solutions  of  radian  heat  transfer 
problems.  Therefore  in  order  to  avoid  such  difficulty  (when  all  such 
b-b  transitions  are  not  fully  resolved)  one  should  assign  a  correction 
term  to  the  photoionization  edge  v  to  yield  an  apparent  edge  v^.  such  as 

VE  =  VE  "  AES/h  (36) 

where  h  is  the  usual  constant  and  AE  is  the  so  called  advance  of  the 
series  limit  (see  Ref.  5  ) ,  and  let  the  photoionization  edge  cross- 
section  ofv^)  for  this  particular  level  be  defined  at  v^..  In  addition 
the  lowering  of  the  ionization  potential  AE^  should  be  included  in 

the  calculation  of  v_  such  as 

E 
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hv. 


-  AE„ 


nic 


(37) 


where  E^,  is  the  ionization  potential  and  E  is  the  excitation  energy 
of  the  level  in  question. 
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DTI HT AKER  and  GAMMA  FUNCTIONS  - 
DEGCRIPTinii  OF  T./O  ALGOL  PROCEDURES 
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Computer  Sciences  R  &  D 

I 

j 

ABSTRACT 
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i 

Thin  report  describes  the  development  and  testing  of  two  computer 
algorithms  Y.'OG  and  GAMMAC  required  in  connection  with  computation  of  spectral 
properties  of  plasma  arcs  being  carried  out  for  Dr,  C.  H.  Church  of  Quantum 
Electronics  R  &  D.  The  algorithm  V/OG  computes  the  irregular  V?hittaker 
function  \J  , l0  for  real  positive  non-integral  q,  and  for  i  a  non-negative 
integer.  GAMMAC  computes  the  function  f*(z)  for  z  real  or  complex.  A  very- 
good  accuracy  has  been  obtained  from  both  of  these  algorithms.  GAMMAC  is  ur  - 
in  DOG,  which  calculates  one  of  the  four  confluent  hypergeometric  function- 
needed  in  the  calculation  of  transitions  between  eigenstates  described  by  th' 

Coulomb  wave  functions. 
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APPENDIX  D 


wHITTAKER  AND  GAMMA  FUNCTIONS  - 
DESCRIPTION  OF  T\»  ALGOL  PROCEDURES 

By  G.  Basi 

Conputer  Sciences  R  &  D 

ABSTRACT 

This  report  describes  the  development  and  testing  of  two  computer 

algorithms  V.'OG  and  GAMMAC  required  in  connection  with  computation  of  spectral 

properties  of  plasma  arcs  being  carried  out  for  Dr,  C.  H.  Church  of  Quantum 

Electronics  R  &  D.  The  algorithm  WOG  computes  the  irregular  VThittaker 

function  U  , /0  for  real  positive  non-integral  and  for  l  a  non-negative 
R  t  <£+-1/  * 

integer.  GAMMAC  computes  the  function  p(z)  for  z  real  or  complex.  A  very 
good  accuracy  has  been  obtained  from  both  of  these  algorithms.  GAMMAC  is  ured 
in  Y.DG,  which  calculates  one  of  the  four  confluent  hypergeometric  functions 
needed  in  the  calculation  of  transitions  between  eigenstates  described  by  th'> 
Coulomb  wave  functions. 
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Summary 


P  r-  des;rlbes  the  development  end  test,_  .  f 

alconthms  :,cc  and  G-Jf  ir  >  ,  *•  t  a  two  computer 

0'*“'C  re‘>ulr'id  t"  connection  with  coroutat, 
properties  of  „  .  Potation  of  spectra 

-  arcs  being  carried  out  for  Lr.  c.  „  rh„  „ 

Electronics  &  D  ^  *  Church  o*  :.uu:  turn 

"  1116  alGorithn  \.'0G  computes  the  t 

function  f,,  ,  the  -regular  Lhittaher 

H.t+yz  r»  real  positive  non-integral  and  e™ 
integer.  ^  a  1  ^n-negutivr- 

”  “  COnPUteC  the  f,JnCtior  r(=)  *  real  or  c 

■“  <«  *en  obtained  from  both  c-  - 

in  WO.  which  calculates  one  of  -be  fo  ^  ^  CW*C  1=  ““ 

needed  in  the  calculate  -  t  ‘  *-*  — 

culation  01  transitions  between  eigen-ate  , 

Coulomb  wave  factions.  °  descrlbed  by  the 

Introduction 

Purpose  of  this  study  was  to  ,  a 

**»  fortran  into  Algol  .„d  lmprove  upcr  "‘J  ”*”**"'  "* 

— -  „d  c^rr:  n'ither  °f  th- 

Theae  algorithms  will  be  uaed  in  a  °"  """  n'C"’S"Jr- 

ed  in  a  program  to  calculate 

the  eig.not.tea  described  bv.br  transition,  between 

y  the  Coulomb  wave  function. 

Ih«  calculation  Of  tb,  transition,  involve,  the  us.  of  confluent 

hypergeometric  functions.  Tuo  of  the  en 


x  H  -  (b-x)  &  .  av  , 
dx^  ux 


^  tW°  SOlUti°nS  °f  —  -  Kummer ' s  functions  F  (..  b.  x)  Md 
(  >  *>;  x)  defined  as  follows: 


y'AUi  biX)'  tWn  £ 

n=0 
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Where  (a)n  =  a(a+l)  -  (a+n-l)  for  n*l  and  (a)  =  l 

o 

(b)n  is  similarly  defined. 


y .  „<»•  b.  x) _ 2 _ j  b;  «> 

Sinrt  i  r(l*a-b)r(b)  '  Wa)r(?-bl - 


The  other  two  confluent  hyper geometric  functions  are  the  solution; 
to  Whittaker's  enuation2 


A,! 


,2  : 

±2  +  .  a  r+i , 

d22  U  +  i  -  —i  y  -  o  • 


z"  / 

Two  solutions  are  Whittaker’s  functions  M^^z)  and  W^  /+l/o(z). 
Whittaker's  functions  are  defined  as  follows: 

Mti,/+1/2(z)  *  zi+1e_1/2z  ^U+l-r^+^z). 

-  n 

i 

WT|,i*i/2(z>  ■  sin(2»i7.-  /FT-  jlfirfgiiy  *  [ 

1.1  Available  programs 


Two  programs  that  calculate  the  generalized  hypergeometric  function 

pFq^al,a2”  ‘ 'ap;  bl'b2' *  * *bq’z)  are  available  in  Algol  B5000. 

These  programs  were  written  by  Hammers3.  One  of  these  two  programs  computes 
a  general  hypergeometric  function  which  has  a  complex  or  real  argument  and 
any  number  of  real  parameters.  The  other  program  computes  the  same  function 
for  any  number  of  real  or  complex  parameters. 
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There  are  three  programs  called  BEAUTY,  BEAST,  and  CHIC  in  Fortran  II. 
Theae  programs  were  written  by  Johnson  and  Sangren4  and  calculate  any  generalized 
hypergeometric  functions.  Program  BEAUTY  can  be  used  to  calculate  any  generalized 
hypergeometric  function  with  up  to  ten  numerator  and  ten  denominator  parameters 
and  one  argument.  Program  BEAST  can  be  used  to  compute  analytical  continuation 
of  an  ordinary  hypergeometric  function.  The  program  called  CHIC  can  be  used 
for  calculating  3k  distinct  hypergeometric  series  in  two  variables. 

Two  other  programs  WOG  and  HLIGAM  are  obtained  from  McGuire5.  Both 
these  programs  were  obtained  in  Fortran.  Program  WOG  calculates  the  irregular 
Whittaker  function  /+1y2(z)  ior  tj>o,  real,  and  non- integral.  Program  HLIGAM 
calculates  the  same  function  where  rj  is  a  pure  imaginary  number. 

Both  of  these  programs  call  for  another  program- -GAMMAC--  which  calcu¬ 
lates  the  function  r(z)  for  any  z  real  or  comples.  GAMMAC  was  obtained  in 
Fortran  from  Dr.  Mcguire. 

2.  Computation  of  GAMMAC 

Procedure  GAMMAC  computes  the  value  of  r(z)  in  the  following  way: 

Let  z  .  (x,y)  ji.e.  z  =  x+iyj.  The  notation  z.  =  (x^)  will  be 
used  during  the  discussion  and  shall  be  understood  for  i  -  1,2, 3,4.  £  wilJ 

denote  the  complex  conjugate  of  any  complex  number  z. 

For  positive  real  integral  values  of  z,  r(z)  is  computed  as  (z-1):, 
for  non-positive  real  integral  values  of  z,  the  procedure  does  not  work. 

If  x<0,  then  the  procedure  uses  the  following  recurrence  relation  to  relate 
T(z)  to  riz^  where  the  real  part  of  zl  is  positive. 

or  r<*>  " 

r(*-y)  '  K.y)  ■ 


v-u 


Now  if  (which  is  the  same  as  y),  in  is  negative, the  rCz^)  is 
related  to  r(z^),  where  y^  is  positive,  by  the  following  relation 

rf^)  -  rdj)  =  r(z2) 

i.e.  z^  =  z^,  or  y^  a  -y^  and  ■  x. 


with  the  use  of  the  above  relations  the  computation  of  r(z)  for  general  z 
reduces  to  the  computation  of  r(z^)  whose  both  real  and  imaginary  parts  are 
positive. 


The  next  step  is  to  relate  r(z^)with  r(zy  where  0  as  y  <  1.  This 
is  done  by  making  use  of  Gauss'  multiplication  theorem.  Gauss'  multiplica¬ 
tion  theorem  asserts  that  if  z„  n  nz,  where  n  is  an  integer  then 


rU2) 


r(*3)r(yl)r(z3*f> . r(z^) 

,p  i  (n-1)  K 

(2  it)  n 


In  our  case  integer  n  is  chosen  to  be  the  next  greater  integer  than  y^  so  that 
y2/n  ■  y  <  1.  Next  we  relate  T(z^)  to  rCz^)  where  z^  is  such  that  both  Xj^ 
and  y^  satisfy  the  inequalities  0  <  x^  i  1,  0  i  <  1.  In  fact  y^  -  y 3  and  x)+ 
-  pc^Jin  case  x^  is  a  nonintegral  real  ~  -od  if  x^  is  an  integer  x^  * 

where  greatest  integer  contained  in  x^.  The  relation  between  T(x^)  and 

r(z^)  can  conveniently  be  expressed  as  follows. 

r(x3)  -  (z.,-i)(z3-2)(z3-3)  .  (z4)r(zu) 

Now  the  computation  of  r(z)  for  any  complex  number  z  has  been  related  to  r(z^) 
where  z^  is  such  that  0  <  i  1,  0  *  y^  <  1.  r(z^)  is  determined  by  t.he 
Padd'  power  approximation  of  l/r(z). 

The  procedure  GAMMAC  was  obtained  in  Fortran  and  has  been  translated 
into  Algol.  Some  changes  were  made  during  the  translation  to  make  the  program 


more  efficient  and  clear.  A  basic  constant  vas  defined  to  greater  accuracy 
to  reduce  the  error.  At  present  it  is  as  discussed  above,  and  calculates 
the  function  r(z)  for  a  real  or  complex.  The  values  Mere  checked  against 

the  available  tables  and  found  to  be  accurate  to  eight  significant  figures 

in  the  region 


x  =  -2.0  (.25)  2.0 

y  =  -2.0  (.25)  -  2.5  and 
y  =  .25  (.25)  2.0 

3‘  £25Putation  of  Irregular  Whittaker  Function  W  (z)  for  r  _eal 

positive  and  non-integral.  ,!kJ+l/p  '  1 - ~~L 


Progi'am  WOG  calculates  the  irregular  Whittaker  function  W  (z) 

r  a  JbJ+l/_2K  ' 

or  tj  real,  positive  a»;d  non-integral.  The  function  is  r(V*+l) 

computed  by  two  different  expressions  for  smaller  and  larger  values  of  2R/t) 

where  3^  -  1,  The  expression  used  to  compute  the  function  for  small  values 
Of  Z  is  given  by  Hartree6  and  is  as  follows. 


,i+l  Jh  i+1/2^  ^  Sin  »  tj  e-.l/2z 


(_)^X  - - 

n  r(n-/) 


(nz)' 


2/ 


r(V/+i-m 


(2^-m)  m 

m  2 


m 


+  S^T  r(T)rii-l)  _ 

rT^f  Z-rl-m)  m*  (m-2i-l)'. 


i  log  Z  +  *  cot  *T)  + 


n  +  *(ri+i+l-m).^(m4-l)-^(m_?/) 


(1) 
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Some  important  aspects  of  the  evaluation  of  this  expression  are  worth 
mentioning.  Computation  of  the  multiplying  factor  makes  use  of  the  identity 


r(Tfr/+i) 


2i+i  r(T)-i) 


A 


(AP) 


In  computing  the  finite  sum  in  Equation  (l),use  is  made 


of  the  identity 


TT 

k-l 


sen 


m 


(A5) 


In  summing  the  infinite  series  in  Equation  (l)  the  summation  parameter  is  trans¬ 
formed  by  k  ■  m  -  2/  -  1  where  k  goes  from  0  to  The  term,  when  k  .  0  is 
obtained  separately.  The  expression  containing  *  functions  is  computed  in  three 

different  ways  depending  upon  the  value  of  tj.  In  case  t,  <  i  and  tj  <  60  the 
identities 

n-1 

f(n)  *  -  y  +  y  I  (n*2) 


k-i 


♦(l+z)  *  -r  + 


z 

k(k+z) 


2  J4  -1»  -2,  -3, 


(2) 


(3) 


^(l-z)  m  l)r(z)  +  *  Cot  JI 


(V 


are  used. 


If  i  <  tj  <  60  the  following  additional  relation  is  used, 
^  ^  **  7n-l)  +  z  +  (n-2)+z  + . +  ^z) 


(5) 


In  the  case  n  >  60  the  asymptotic  expansion  of  the  *  function  is  us 


ed, 
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♦  (z)~ln(z) 


22  12 z2 


120  z  252z 


(z  “  in  |  arg  z  |  <  r) . 


When  the  argument  z  of  the  Whittaker  function  is  large,  the  asymp¬ 
totic  expansion  is  used  for  the  computation  in  place  of  Equation  (l) .  This 

Q 

expansion  is  given  by  Whittaker  and  Watson  . 


S,/+l/2 


"I2 

(z)~e  z 


00 

Hal 


i+l/2)  -(n-1/2) 


n!  z 


/+1/2) 


Some  test  runs  with  McGuire's  program  gave  results  of  rather  low 
accuracy- -two  or  three  significant  figures,  and  some  steps  were  taken  to  im¬ 
prove  this.  The  major  source  of  error  appeared  to  be  the  summation  for  the  * 
function  in  Equation  (3),  which  is  slowly  convergent.  McGuire  truncated  this 
summation  at  50  terms.  We  have  summed  greater  numbers  of  terms  and  also  intro¬ 
duced  a  correction  for  the  remainder  term  by  consideration  of  the  integral 


oo 

f  rV  -  i  W-SU 

/  x(  Xa)  a  |xo-a } 


By  approximating  to  the  integral  in  two  different  ways,  we  determine  the  errj' 


bounds  for  the  remainder 


1.  N-l 
a111  N-l-a 


V  1  '  1  In  N'l/g  \ 

L  "  N-l/ 2 -a  j 

J-N 


The  remainder  term  is  approximated  by  the  mean  of  the  two  bounds  and  good 


accuracy  is  obtained  by  taking  N  -  1001.  In  the  case  of  t;ie  infinite  oun  involvea 
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m  (1),  the  series  .as  belnS  truncated  .hen  the  tern  became  less  than  if,'3 

of  the  sum.  in  the  .edified  version,  the  suction  is  truncated  when  the 

tern  becomes  s, nailer  then  lo^  of  the  „  and  .hen  at  least  10  tents  have 
been  summed. 

m  computing  the  asymptotic  expansion  of  Whittaker's  function,  the 
number  of  terms  summed  by  McGuire  .as  the  greatest  integer  contained  m  ,  +  , 
The  reason  for  this  choice  is  not  clear  and  It  .as  found  tut  improvement  could 
he  obtained  by  string  up  to,  but  excluding,  t„e  smallest  term  (norwl  procedur 
with  en  asymptotic  series).  This  method  us,  therefore,  adopted. 

With  these  modifications,  a  considerable  increase  in  accuracy  has 

teen  obtained,  at  the  expense  of  some  increase  in  computer  time.  Irbies  of  the 

Whittaker's  Irregular  function  vere  not  available.  Por  checking  the  accuracy 

of  the  results  Coulomb  tables^  vere  used  rvmi -mV,  +  vi 

USed*  Coul3rab  tables  list  the  functions 

■*(-£,**)  ^  ^  relation  between  Whittaker's  function  and  the 

functions  P^R)  and  Q{ipR)  is  glVBn  by  the  equation 


w 


(2R) 


(2i+i): 


Q^^R)  Sin«i]  +  p  (i  ,r)  Cosjftj 


results  of  the  prog™  have  been  found  to  a^ee  with  the  results  obtained 
ty  using  the  tabulated  values  to  six  or  more  significant  figures  for  the 
ranges  /  -  0,  1,  2,  l/,2  „  -2.0(.2)-.2,  ,  .  .5(.5)k.O.  ftom  result,  obtalDed 

in  these  ranges,  it  appeared  appropriate  to  switch  the  computation  from  the 
series  (1)  to  the  asymptotic  form  (7)  at  the  value  t  „  9. 
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th.  cal  uxT  ““  St'P  ^  th'  W0JeCt  ^  Urltl°8  "  C<*PUl"r  Pr0era” 
thp  calculation  of  trans-f  *  in„e  u  * 

raault  ,  ™  «»  to  check  the 

-or  th.  o.her  Progrsa  called  HL1GAM  against  the  available  Coulanb 

Z Z  ““  "  P0SSlWe  *—  lf  *  —  Its  accuracy, 

e  ure  HLIGAM  calculates  the  function 

(2r\ 

ru+l+itj)  WiT),i+i/2  ( Irf  / 

^  °"ber-  ^  <«. 

P  obliti  ^  “  rtU  "  iD“^  *  «.  transitions 

9  Integration  vlil  be  checked  against  the  tables 

given  by  Peachy. 
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Program  Usage 


Mx  GAMMAC  (x,y,u,v)j  An  aigol  procedure  to  compute  r(z)  for  real  or  complex 
values  of  the  argument  z. 

Procedure  GAMMAC  (x,y,u,v)  computes  the  function  !(z)  for  any 
argument  2,  real  or  complex.  This  procedure  can  be  included  in  any  main 
program  by  including  the  following  two  cards  in  the  progrume  declarations 

$$  A  GAMMAC 

COMMENT  HIGHER  SEQUENCE  CARD  ; 

Both  the  cards  must  have  some  sequence  number  In  the  column  73-80.  and  the 
sequence  number  on  the  comment  card  m-t  be  higher  than  that  on  the  $$  ft  GAMMAC 
card.  The  fontal  parameters  x  and  y  are  the  input  parameters,  x  and  y  are 
b.tb  real  and  re-pectively  are  the  real  and  imaginary  part  of  the  complex 
number  t,  whose  gamma  function  Is  required,  u  and  v  are  the  output  parameters 
and  respectively  are  the  real  md  imaginary  parts  of  the  complex  number  r(z). 

For  example  if  one  wants  to  calculate  r(t.)  and  store  the  value  at 
z2  vtere  =-  (xj.y^  for  1  -  1,2,  are  complex  numbers  then  one  would  do  so 
by  including  a  card,  at  the  desired  place,  in  the  main  program  with  the 
following  information  on  it. 

GAMMAC  (x1,y1,x2,y2); 

Either  one  or  both  of  the  input  parameters  x  and  y  can  he  real  valued  arithmetic 
expressions.  There  Is  no  loss  of  efficiency  caused  by  using  the  expression  as  an 
acVrnl  input  parameter  since  it  is  called  by  value. 
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Limitations.  The  procedure  has  the  following  limitations. 


For  the  non-positive  integral  values  of  the  argument  an  error 
message  of  divide  by  zero  would  occur  since  the  gamma  function  is  singular 
at  such  points.  For  the  rest  of  the  negative  real  values  of  the  argument 
procedure  calculates  I'(z)  except  very  close  to  the  negative  integers  where  an 
exponent  overflow  message  would  occur.  For  positive  real  values  of  the 
argument  procedure  computes  1(2)  for  z  53-3-  For  larger  values  of  z 
an  error  message  of  exponent  overflow  occurs.  For  very  small  values  of  z 
(z  <  8.758  @-47)  an  error  message  of  divide  by  zero  occurs  since  computer 
sets  such  numbers  to  zero. 

For  pure  imaginary  argument  |z|  >  66. w  the  value  of  |r(z)|  is 
so  small  that  it  is  set  to  zero  by  the  computer. 

For  z  =  a  complex  number  on  the  digonal  x  *  y  or  x  *  -y  and  x  > 0. 
While  computing  f(z)  of  z  such  that  >  >  42.00  an  error  message  of  maxinum 
argument  of  exponent  would  occur.  On  the  diagonals  x  =  y  and  x  =  -y  for 
x  <  0  the  value  of  r(z)  is  so  small  for  x  <  -22.60  that  it  is  set  to  zero  by 
the  computer. 
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WOG  (E,  ETA,  R,  W)j  An  algol  procedure  to  compute  the  irregular 

V 

Whit  taler  function  n,fc-H/2v  r,; 

r{ n+it+i) 


The  procedure  WOG(L,ETA,R,W)  calculates  the  function 


w  (ii) 

n.l+l/2(  n 

r(n+e+i) 


A  deck  of  «-.:.rdE  for  this  procedure,  tc  be  included  in  the  main  program, 
may  be  obtained  from  the  writer.  The  formal  parameters  L,  ETTA  and  R 

correspond  to  the  variables  n  and  R  in  the  Whittaker  function  w  .  (— ) 

n,*+l/2  n 

Either  one  or  any  number  of  these  input  variables  L,  ETA,  and  R  can  be 
a  real  valued  arithmetic  expression.  There  it  nc  loss  of  efficiency  caused 
by  using  an  expression  as  an  actual  parameter  since  it  is  called  by  value. 

W  is  the  output  parameter,  which  has  the  value  of  the  function 
as  computed  by  the  procedure. 

One  of  the  infinite  series  involved  in  computation  is  a  function 
of  Z  and  ri  aut  not  R.  This  summation  is  the  most  time-consuming  part  of  the 
computation.  To  increase  efficiency  when  successive  values  of  W  are  to  be 
calculated  with  constant  1  and  n  ,  varying  R,  the  procedure  has  been  written 
so  that  this  particular  sum  is  preserved  between  calls.  It  is  recomputed  only 
when  either  torn  changes.  Therefore  it  is  reconmended  that  the  user  arrange 
his  program  so  that  consecutive  calls  on  the  procedure  WOG  are  made,  as  far 
as  possible,  with  constant  l  and  n  (i.e.  variation  of  R  should  be  the  inner¬ 
most  loop) . 

Limitations.  This  procedure  is  limited  to  the  real  values  of  R,  to  n 
positive  and  nonintegral,  and  to  l  integral  only. 
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APPENDIX  E 


RADIATIVE  TRANSPORT  IN  A  XENON  ARC 
by 

B.  W.  Swanson 


I.  INTRODUCTION 

The  problem  of  interest  has  been  to  determine  the  steady  state 
temperature  profile  in  a  Xenon  arc,  which  requires  a  solution  of  the 
energy  equation 


9  S  1  3S  ■+  ■*  2 

_  +  v  .  F  +  0(PjS)E;  ,  0 


where  F(r)  is  the  radiative  flux  vector  T  M 

uj..  \ector,  i  che  plasma  temperature,  o  the 

electrical  conductivity,  P  the  pressure  and  E  the  electric  field.  The 

Schmitz  function  S  is  defined  by  the  equation 

T 

S(T)  =  J  K(T)dT  (2) 

TW 


and  the  divergence  of  F  is  given  by 


V  • 


(3) 


where  the  radiative  flux  F  is  found  from  the  equation 

00 

J  1^  (r  ,(i)lu)dwd  a 

o  ui*47t 


(4) 


In  equation  4,„  is  ,  unit  direction  vector  and  Ij(r,S  Is  the  monochronu 
tic  intensity  of  radiation. 


E-l 


The  first  phase  of  this  program  consisted  of  writing  a  compu¬ 
ter  program  for  calculating  the  radiative  flux*  which  can  be  expressed 

as 

Amax  it 

'“’■"*// 1  K^(s)  B^(s)  G^tg^sJJcos  <p  dsdifidA  (5) 

Amin  o  o 


In  equation  5,  K^(T)  is  the  spectral  absorptivity,  Amin  and  Amax  the 
wavelength  band  of  interest,  B  is  the  Planck  function, 


and  the  function  G  is  given  by 

JL  x 

2  -  sin0 

V’  ' 


(sine)  d6 


Kx(c)d; 


(6) 


and  accounts  for  the  attenuation  of  radiation  by  self  absorption.  The 
geometry  for  the  evaluation  of  equation  S  is  shown  in  figure  1. 

Since  the  numerical  evaluation  of  the  triple  flux  integral  is 
expensive  the  second  phase  of  the  program  sought  a  faster  flux  approxi¬ 
mation  method.  An  analysis  was  made  which  permitted  an  a'priori  integra¬ 
tion  with  respect  to  wavelength  to  reduce  computation  time.  The  approx- 

2 

imate  flux  integral  is  given  by  the  equation 


F(r) 


V*’ 


-*  J  J  Vs> B(s)  Gi  1 


■V»> 


K  (Odc]cos  0  dsd$ 

a 


(7) 


where 


B(s) 


B,(T(s)]dA 


(8) 


and  the  mean  absorption  coefficient  K  (s)  is  defined  by  the  equation 


E-2 


(9) 


T  (s ) 

Vs>  "  'FTT^T'  Vs>  +  V8> 


The  term  b  in  a  constant  of  the  order  of  unity  which  can  be 

varied  to  improve  the  approximation.  The  terms  t  and  t_  are  the 

R  P 

Rosseland  and  Planck  optical  lengths  along  the  ray  i'rom  point  M  to 
point  R  in  figure  1  and  are  given  by 

s 

-  I  ^(OdC 

0  (10) 

■r 


.R(s) 


Tp(s) 


Kp(OdC 


where  the  Planck  and  Rosseland  absorption  coefficients  are  given  by 


Ka(T)  Bx(T)dA 
B(T) 

1 

K  (T)  dT  ’  dA 

A 

dB  ” 
dT 


(ID 


and 


B(T) 


2 

Bx(T)dA 


(12) 


Tbs  form  of  the  approximation  is  such  that  under  optically 
thin  conditions  when  td  and  t  are  much  less  than  unity,  the  mean  absorp- 
tion  coefficient  approaches  the  Planck  absorptivity  Kp  and  equation  7 
is  exact.  Under  optically  thick  conditions,  when  td  and  t  are  much 
greater  than  unity,  approaches  the  Rosseland  mean  and  equation  7 
is  a  good  approximation. 
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The  subject  of  this  report,  which  is  the  third  phase  of  the 
program,  consists  of  the  numerical  solution  of  equation  1  utilizing  the 
exact  and  approximate  flux  programs. 


II.  XENON  PROPERTIES 

The  initial  data  needed  for  a  numerical  solution  consists  of 
the  electrical  conductivity,  thermal  conductivity  and  spectral  absorp¬ 
tivity  as  functions  of  temperature  and  pressure.  The  temperature  de¬ 
pendence  of  the  S  integral  is  shown  in  figure  2  for  a  pressure  of  11 
atmospheres,  and  Lhe  electrical  conductivity  o  is  shown  as  a  function 
of  S  in  figure  3.  The  variation  of  spectral  absorptivity  with  tempera¬ 
ture  is  shown  in  Table  I.  An  inspection  of  the  table  reveals  that  the 
spectral  absorptivity  in  the  ultra  violet  region  is  orders  of  magnitude 
greater  than  that  in  the  visible  region.  Therefore,  the  flux  computer 
program  consists  of  an  exact  flux  program  for  the  U.  V.  and  an  approxi¬ 
mate  flux  program  for  the  visible.  The  wavelength  bands  defining  these 
regions  are  given  by  the  equations 


U.V. 


Visible 


Amin  ■  .0612  X  10 
max  *  .0952  X  10 


.1034  X  10 
2.0  X  10~*  cm 


cm 

cm 

cm 


III.  NUMERICAL  ANALYSIS 

A  steady  state  solution  of  the  energy  equation  can  be  obtained 
by  solving  the  "transient'’  equation 


3_S  _afs  1  3S 
3t  "  dr2  r  3r 


ldf  +  f)  +  °(s>e2 


(13) 


until  a  steady  state  solution  is  achieved.  The  following  calculation 
procedure  was  followed: 
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\ 


Table  I  - 

Spectral  Absorptivity  of 

Xenon  v 

-K 

:m  ) 

4 

AxlO 

cm 

i 

T 

l°Kl 

.  3000 

6000 

>  10000 

140C0 

1  18000 

2.0 

1  1.257  3-14 

4.9090-05 

.2502 

3.046 

2,406 

1 

t 

1 . 7500 

I  2.6160-14 

1.6790-05 

.0523 

.5612 

.4115 

i 

i 

.4615 

.  1.2240-12 

1.0520-04 

.1352 

.9441 

.5476 

i 

• 

) 

.3333 

8.184  0-12 

.0768 

.5461 

.3027 

« 

i 

i 

.2609 

6.7750-11 

1.076@-04 

.0511 

.3594 

.1872 

1 

t 

tu 

.2143 

5.5570-10 

1.1940-04 

.0345 

.2446 

.1193 

1 

1 

t 

CQ 

M 

.1818 

9.8420-10 

9.5060-05 

.0225 

.1598 

.0761 

“ 1 

j 

on 

i 

> 

.1579 

6.4460-10 

6.2260-05 

.0147 

.1049 

.0501 

.1395 

4. 449r-J0 

4.2970-05 

.0107 

.0725 

.0347 

.1250 

3.1980-10 

.  3.0890-05 

.00731 

.0521 

.0250 

i 

.1132 

2.3760-10 

2.2950-05 

.00543 

.0387 

.0186 

■  ' 

1.8130-10 

.00414 

.0295 

.0142 

,0952 

2.269.0 

1.134.0 

614.1 

145.9 

13.07 

.0882 

2.094.0 

1.047.0 

56.  .0 

135.7 

12.07 

H 

W 

.0822 

1.920.0 

959.7 

519.8 

124.4 

11.07 

o 

hH 

> 

.0769 

1.746.0 

872.6 

472.6 

113.1 

10.06 

£ 

.0723 

1.572.0 

785.5 

425.4 

101.8 

9.06 

2-1 

nJ 

O 

.0682 

1.397.0 

698.4 

378.3 

90.51 

8.05 

.0645 

1.223.0 

611.3 

331.1 

79.22 

7.05 

.0612  1.049.0  524.2  183.9  67.94  6.04 


(a)  specify  arc  radius  and  central  core  temperature, 

(b)  specify  an  initial  temperature  profile 

(c)  calculate  a  radial  flux  profile 

(d)  define  E  by  the  equation 


o(o) 


(e)  select  a  time  step  At  and  solve  for  a  new  S  profile 

(f)  repeat  steps  c,d,  and  e  until  the  solution  converges 
to  a  steady  state 

The  field  E  is  so  defined  as  to  keep  the  central  temperature 
essentially  constant.  This  procedure  was  adopted  to  avoid  large  central 
temperature  variations  which  would  greatly  affect  the  flux  and  its 

divergence.  The  fiuite  difference  equation  used  to  replace  equation 
13  is 


«*  ♦  ♦ « - 


.  2S”+1] 


*  1 - hir - +  +  ”(SJ)(E”) 


where  the  subscript  i  denotes  a  radial  position.  Equation  14  is  an  im¬ 
plicit  difference  equation,  which  for  linear  equations,  is  known”*  to  be 
unconditionally  stable.  Now  let  0  -  and  define  the  coefficients 
j  by  the  equations 


'i,i-l  0(1 
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Then  equation  14  can  be  written  as 


c  cn"^  +  r  cn+1  =  Rn 

Li,i-1  Si-1  Ci,i  Si  Ci,i+1  1 


(16; 


where 


.  2  F”  .  -  F®  F? 
d n  __  /  -*  n  \  ✓ ,, »  r  ,  i  -1 

Bi£°(ai)(E'  -[  - +  —  ] 


(17) 


Let  S°  and  Bn  define  the  vectors 


form  as 


ft  \ 


r 


n  \ 

2  2  1 


r  =  (  ? 


(18) 


'  sn 

\  N 


Then  the  finite  difference  equations  can  be  written  in  matrix 


+TYV1 

<C>  <>  * 


(19) 


The  row  of  <C>  contains  only  the  oefficients  C  „  . 

K  9 K— 1 9 

C  „  C  and  the  matrix  is  tri-  diagonal.  The  matrix  <C>  can  be 

"inverted"  very  efficiently  by  Gaus's  elimination  method.  The  first 
equation  (i^l)  is  used  to  eliminate  from  the  second  equation  (i*!2), 

the  new  second  equation  used  to  eliminate  from  the  third  equation 

and  so  on  until  finally  the  new  last  but  one  elation  can  be  used  to 
eliminate  from  the  last  equation,  giving  only  one  equation  with 

only  one  unknown  ~.  The  unknowns  S^+E,  S^+E, . . . can  then  be  found 
in  turn  by  back  substitution.  The  method  provides  a  fast  computational 
procedure  for  large  numbers  of  equations.  In  the  numerical  solution. 
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50  raaial  increments  were  employed.  The  flux  program  was  used  to  eval¬ 
uate  the  flux  at  9  points  and  intermediate  values  were  obtained  with  the 
use  of  a  quadratic  interpolation  program. 


IV.  NUMERICAL  SOLUTIONS 

Numerical  solutions  were  obtained  for  n  pressure  of  11  ..tm. 
and  central  arc  temperatures  of  10,000°K  and  18,000°K.  For  the  1 irst 
example,  several  iterative  S  profiles  are  shown  in  figure  4.  Sq  denotes 
the  initial  profile,  and  S^  denotes  the  final  profile  obtained  after  29 
time  increments.  The  final  solution  has  a  maximum  error  c  of  .03%  where 
e  is  defined  as 

15 

£  r  ^2  (20) 

oE 


In  figure  5  are  shown  several  iterative  flux  profiles.  For 
the  final  S  profile  (S^).  the  approximate  flux  calculations  were  com¬ 
pared  with  exact  calculations  over  the  range  anc*  were  found  to 

agree  to  within  1%.  In  figure  6  is  shown  the  final  temperature  profile 
corresponding  to  02^.  The  temperature  profile  is  very  flat  with  a 
narrow  conduction  laysr  near  the  tube  wall.  The  total  current  was  evalu¬ 
ated  from  Ohm’s  law  as 

rR 

I  -  2ttE  /  ro[S(r)]dr  (21) 

“'o 

In  figure  7  is  given  the  variation  of  E  with  Lme",  showing  a  conver¬ 
gence  to  20.69  volts/cm. 

The  iterative  S  profiles  for  the  second  example  are  shown  in 

figure  8.  A  total  of  80  time  increments  wa3  required  for  convergence 

and  the  final  profile  SDn  differs  markedly  from  the  initial  S  profile. 

ou  o 

The  final  solution  Sgg  has  a  maximum  error  of  0.2%.  In  figure  9  are  shown 
the  initial  and  final  flux  profiles.  A  comparison  of  approximate  and 
exact  flux  calculations  over  the  range  ^^“^2  8ave  agreement  to  within 
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must  still  be  considered  an 


10%.  Therefore  the  "final"  solution  Sg^ 
approximation  to  the  exact  solution.  The  temperature  profile  corre¬ 
sponding  to  Sg0  is  shown  in  figure  10.  Lowke  and  Capriotti5  have  shown 
that  the  curvature  of  the  temperature  in  the  central  core  region  is  due 
to  self  absorption  effects.  Also,  the  thermal  conduction  edge  is 
smaller  than  in  the  first  example.  The  iterative  variation  of  E  is 
given  in  figure  11  showing  convergence  to  a  final  value  of  68.95  volts/cm 
and  the  corresponding  current  is  5795  amps.  In  summary,  these  examples 
indicate  that  the  numerical  method  employed  is  stable  and  provides  a 
convenient  method  for  obtaining  steady  state  solutions  to  the  energy 
equation.  Furthermore,  with  slight  modification,  the  method  can  be  used 
to  calculate  actual  thermal  transients  in  radiative  plasmas,  for 
applications  where  transient  effects  are  important. 
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APPENDIX  F 


Analysis  of  Xenon  Discharges 


J .  J .  Lowke 

Westinghouse  Research  Laboratories 
Pittsburgh,  Pennsylvania  15235 


I.  INTRODUCTION 

Calculations  have  been  performed  deriving  temperature  profiles 
of  discharges  in  xenon  flash  tubes  from  given  miterial  functions  for  xenon, 
at  a  specified  experimental  value  of  current  and  electric  field  strength. 

The  calculations  include  rhe  effects  of  conduction,  radiation,  and  self 
absorption  of  radiation.  Comparisons  are  made  of  three  different  methods 
of  doing  these  calculations.  Results  of  these  calculations  indicate  that 
for  typical  values  of  electric  field  strength  and  current  the  temperature 
profiles  are  radiation  dominated  if  the  pressure  is  equal  to  or  higher  than 
four  atmospheres,  the  profiles  are  relatively  fiat,  the  influence  of 

O 

radiation  at  wavelengths  less  than  1000A  is  negligible,  and  the  effect  of 
self  absorption  in  the  visible  continuum  is  of  only  secondary  importance  in 
determining  the  central  temperature. 

An  analysis  has  been  made  of  the  time  dependence  of  the  experimental 
data  of  electric  field  strength,  current  and  radiation.  The  radiation 
measurements  suggest  that  the  xenon  discharge  should  be  treated  as  a  ti 
dependent  problem  rather  than  a  wall  stabilized,  steady  state  which  is 

the  basic  of  the  present  calculations. 
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II.  TEMPERATURE  PROF  .LE  CALCULATIONS 

Calculations  of  the  temperature  profile  for  a  wall  stabilized 
steady  state  arc  for  a  gi  -*n  -iectric  field  strength,  current  and  gas 
pressure  have  been  made  using  (1)  the  Eddington  Approximation,  (2)  a  pseudo 
time  dependent  equation  and  an  exact  radiation  flux  calculation  and  (3)  a 
pseudo  time  dependent  equation  and  an  exact  radiation  intensity  calculation. 

It  is  considered  that  method  (3)  is  the  least  expensive.  However,  for 
comparison  of  calculated  values  of  central  temperature  and  arc  radius  with 
experiment,  a  simpie  constant  temperature  model  is  adequate. 

The  jt  material  functions  for  the  calculations  are  those  supplied 
oy  C.  H.  Church.  The  speccral  absi . btivicies  are  an  estimate  of  the  continuum 
and  neglect  line  radiation.  Thermal  and  elect  a  conductivities  are  calcu¬ 
lated  using  a  program  of  De  Voto.  The  effects  of  radiation  were  determined 
by  dividing  the  spectrum  into  20  equal  frequency  bands  extending  to 
5  x  1016Hz(600A). 

1.  Method  using  Eddington  Approximation 

This  method  is  outlined  in  reference  (1).  The  following  equations 

are  solved  numerically. 


\ 
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V. (F  +  Fd)  -  oE 
C  K 


(i) 


Fc  »  -kVT 


(2) 


Fd  dv 
Rv 


(3) 


V.Fd  -  4tt  (£  -  J  K  ) 

Rv  v  v  v 


-►  4r. 

FRv  VJv 


(A) 

(5) 


The  nomenclature  is  given  in  Table  I. 

Equation  (5)  is  the  only  inexact  aquation  and  is  known  in  astro¬ 
physics  as  the  Eddington  approximation.  The  equation  Is  in 
effect  a  diffusion  approximation  for  the  radiation  flux  and  is  equal 
to  a  constant  multiplied  by  the  gradient  of  the  radiation  energy 
density,  which  is  proportional  to  VJ^  .  An  iterative  procedure  is 
used  to  satisfy  the  boundary  conditions  which  are  as  follows 


at  r  »  0 


F  -  0 
Rv 


snd  F_  ■  0 
X 


T  ■  T, 


W 


end  j  *•  — r  F_ 
v  2  Rv 

TT 


at  r  =  R 


Table  I 


Nomenclature 


Conduction  Flux  Density 

Radiation  Flux  Density 

Thermal  Conductivity 

Electrical  Conductivity- 

Spectral  Absorbtivity 

Spectral  Emission  Coefficient 

Radiation  Flux  Density  per  unit  Frequency 

Radiation  Intensity  per  unit  Frequency 

Temperature 

Electric  Field  Strength 
Electric  Current 
Radius 

Wall  Temperature 

Specific,  lleat  at  Constant  Pressure 
Density 

Integrated  Emission  Coefficient 
Time 

Planck  Function 


—2 

(watts  cm  ) 

_o 

(watts  cm  ) 

(watts  cm_1°K_1) 

(mho  cm  *) 

(cm  1) 

-3  -1 

(watts  cm  ster  sec) 

(watts  cm  ‘'sec) 

-2  -1 

(watts  cm  ster  sec) 
( °K) 

(volt  cm  ) 

(amps) 

(cm) 

CK) 

(Joules  gm  ^°K  ^) 

(gm  cm  J) 

(watts  cm  Jster  1 ) 
(sec) 

-2  -1 

(watts  cm  ster  sec) 


Pressure 


(Attrospheres) 
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Input  values  are  the  material  functions  k,  a*  and  e^  (which 
are  functions  of  temperature),  end  the  experimental  values  of 
I  and  E.  R  is  determined  from 

2rR(Fc  +  FR)r=R  =  IE  (6) 

The  solution  of  these  equations  giving  T,  FR  and  Fc  as 
a  fuaction  of  radius  for  E  *  29  volt/cm  and  I  *  1140  amps  at  a 
pressure  of  11  atmospheres  is  shown  in  Fig.  1. 

2.  Method  using  a  pseudo  time  dependent  equation  and  an  exact 
radiation  flux  calculation 

This  method  has  been  suggested  by  B.  Swanson. 

Rather  than  iterate  to  sat: jfy  split  boundary  conditions 
as  in  the  former  method,  an  assumed  initial  temperature  profile 
is  modified  using  che  time  dependent  equation  to  obtain  the 
steady  state  solution. 

The  time  dependent  equation,  neglecting  effects  of  mass 
flow  is  given  by  equation  (7) 

3T  2  -*■  -* 

CpP-  =  oE  -V.FC-  V.Fr 


(7) 
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As  we  are  only  interested  in  obtaining  the  sterdy  state  solution, 

C  p  is  made  equal  to  1,  a  constant.  Fr  is  obtained  using 
P  k 

—V 

equation  (2)  and  F^  is  obtained  exactly  using  the  integration 

2 

procedure  of  Swanson.  A  number  of  difficulties  associated  with 

this  method  are  overcome  as  described  below. 

(1)  The  numerical  solution  of  this  equation  is  unstable 

for  large  step  s^.zes  in  time.  This  instability  is  described  by 
3 

Rich^myer  and  in  the  present  method  of  solving  the  equation,  the 
size  of  the  time  step  for  stability  has  been  empirical]’'  found 
to  closely  follor  the  requirement  cited  by  Richtmyer  that 


where  F  *  k/pC 


Generally  the  small  step  size  requirement  makes  computing  times 

prohibitively  large.  However  in  radiation  dominated  arcs,  the 

temperature  of  most  of  the  arc  is  controlled  by  the  terms 
2  * 

oE  -  '*F  of  equation  (7)  and  not  by  VF  which  causes  the 

K  C 

instability.  By  Initially  dividing  the  termV.F^  by  a  factor  of 
the  order  of  30  and  increasing  the  step  size  by  the  same  factor 
the  rate  of  convergence  to  the  final  profile  is  greatly  increased- 
The  multiplying  factor  is  then  removed  to  determine  the  outer  small 
region  of  the  profile,  which  is  controlled  by  conduction.  In  this 
latter  stage  of  the  calculation  it  is  only  necessary  to  recalculate 
F^,  whi"h  is  the  most  time  consuming  part  of  the  calculation,  at 
say  every  twentieth  time  s^ep,  for  Fn  is  determined  primarily  by 
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the  temperature  of  the  inner  region  of  the  arc.  With  this 
simplification  any  initial  temperature  profile  can  be  made  to 
rapidly  converge  to  the  final  solution. 

(2)  In  the  present  problem,  current  and  electric  field 
strength  are  input  values  and  the  radius  and  the  temperature 
profile  are  the  output.  After  each  time  step  the  radius  is 
recalculated  using  equation  (9) ,  where  the  temperature  profile 
is  regarded  as  a  function  of  r/R 


R 


2 


I  [2irE 


1 

xodx] 


j 

o 


(9) 


where  x  =  r/R. 

In  this  way  the  correct  current  is  maintained  as  the 

temperature  profile  is  varied, 

(3)  It  is  found  that  the  temperature  profile  can  vary 

monotonical ly  and  very  slowly  in  a  way  making  it  difficult  to 

determine  the  error  from  the  true  steady  state  solution.  To 

test  for  convergence,  after  calculating  each  new  temperature 

profile,  equations  (1)  and  (2)  are  solved  to  obtain  F  as  a 

R 

function  of  r.  If  the  temperature  profile  is  correct,  the  values 

■>  •¥ 

of  F  will  agree  with  the  values  of  F  obtained  from  the  spectral 
R  R 

absorbtivities  using  the  Swanson  exact  integration  procedure. 
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In  the  numerical  results  presented  using  this  method  the  two 
values  of  FR  agree  to  U  in  each  case.  The  method  was  used  for 
the  conditions  of  Fig.  1  and  the  temperature  profile  is  shown 
as  the  broken  curve.  The  small  difference  in  the  profiles  near 
the  outer  wall  is  considered  to  be  due  to  an  insufficiently  small 
step  size  in  method  2  and  not  due  to  the  Eddington  approximation. 

The  step  size  in  r/R  near  the  wall  in  method  2  was  0.01  whereas 
in  the  Eddington  method  it  was  0.0001. 

The  method  was  also  used  for  p  =  4  Atmospheres  with 
E  =  29  volt  cm  and  I  =  1140  amps  giving  results  insignificantly 
different  from  those  shown  in  Fig.  2,  which  were  obtained  using 
method  3.  Fig  3  gives  results  obtained  using  method  2  showing 
the  variation  of  FRy  with  r  at  two  different  frequencies.  The 
radiation  flux  integrated  from  the  absorption  edge  at  v  =  2.9  x  10i5 
to  v  =  5  x  1015  Hz  attains  a  maximum  of  only  2%  of  the  total 
radiation  flux  at  r/R  =  0.87  for  the  conditions  of  Fig.  2. 

The  inclusion  of  lines  in  this  region  of  the  spectrum  will  only 
decrease  the  radiation  flux  as,  from  the  Eddington  Approximation, 
which  is  equation  (5),  Ky  is  increased  in  the  region  of  lines  and 
Jv  '  Bv  (see  FiS-  consequently  it  is  concluded  that,  except 

for  very  h  gh  powers  where  the  central  temperature  is  greater  than 
14000  K,  the  effect  of  radiation  for  wavelengths  less  than  1000A 
can  generally  be  neglected  in  the 


present  calculations. 
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3.  Method  using  a  pseudo  time  dependent  equation  and  exact 
radiation  intensity  calculations 

For  the  conditions  of  the  present  problem  it  is  of 
considerable  advantage  to  express  equation  (7)  in  the  form 


C 

P 


P 


9T 

dt 


oE2  -  V.F  -  4ir(e 
c 


J  K  dv) 

V  V 


(10) 


where  use  is  made  of  the  continuity  equation  for  radiation, 


B  K  dv  and  c  is  simply 

v  v 


equation  (4).  From  Kirchoff's  Law  e  ■ 

a  function  of  temperature.  can  be  obtained  exactly  from  a 
modification  of  the  Swanson  integration  procedure  used  to  obtain 


The  advantages  of  this  modification  to  method  2  are 

a.  The  increment  in  temperature  does  not  depend  on 

the  gradient  of  a  computed  quantity.  In 

method  2  to  determine  7.FD,  F  needs  to  be 

R  R 

calculated  to  high  accuracy  and  at  many 
-adial  steps. 

b.  For  typical  experimental  conditions  in  xenon, 

J^K^dv  is  small  compared  with  e,  provided  that 
ultraviolet  radiation  can  be  neglected.  Then 
the  integral  need  only  be  calculated  at  say 
every  100th  step,  instead  of  after  every  20th 
step,  as  outlined  in  the  latter  part  of  Section 
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II  2  (1).  In  fact  if  the  integral  in  the. 
visible  range  of  the  spectr  ,m  is  set  at  zero, 
the  central  temperature  is  decreased  by  only 
200°K  for  Fig.  1  and  300°K  for  Fig.  2. 

Instead  of  having  E  and  I  as  the  input,  as  in  method  2, 

E  and  R  are  treated  as  the  input  and  the  value  of  I  is  calculated 
as  output.  The  difficulty  outlined  in  II  2  (2)  is  thus  removed 
and  the  stability  of  the  numerical  solution  is  also  improved. 

The  profile  obtained  using  this  method,  with  E  ■  29 
and  R  ■  0.  1  as  input  at  a  pressure  of  4  Atmospheres  is  shown 
in  Fig.  2.  Fig.  4  gives  as  a  function  of  frequency  at  two 
radial  points. 

4.  Constant  temperature  model 

In  any  comparison  of  theoretically  derived  temperature 
profiles  with  experiment  one  of  the  principal  interests  is  to 
compare,  for  any  value  of  I  and  E,  the  theoretical  central  tempera¬ 
ture  and  radius  with  experimental  observations.  An  assessment 
can  then  be  made  of  the  accuracy  of  the  input  material  functions,  c 
and  ,  assuming  that  the  experimental  arc  is  wall  stabilized  and 
in  an  effective  steady  state  condition.  For  this  purpose,  for  the 
experimental  conditions  considered  in  this  paper,  the  central 
temperature  and  R  can  be  computed  to  with.'.  5%  simply  by  assuming 
that  the  arc  is  at  constant  temperature,  without  computing  the 
temperature  profile. 


II 


ior  radiation  dominated  arcs,  loss  of  energy  by  conduc¬ 
tion  at  the  center  of  the  arc  is  negligible  compared  with  radiation 
losses.  This  the  innut  electrical  energy  equals  the  net  emission 
of  radiation  energy,  as  expressed  in  equation  (11),  i.e., 


aE 


2 


4tt(c 


f 

i  JvKvdv) 


(11) 


At  the  center  of  a  cylindrical  plasma  of  radius  R  and  uniform 
temperature  T,  Jy  can  be  calculated  using  the  two  radiation 
equations  (4)  and  (5)  again  given  below 


7,Fpv  "r  4tt(c  -  J  K  ) 

KV  V  V  V 


4tt 

3K 

v 


VJ 

v 


As  and  Ky  are  constant  the  equations  can  be  solved  analytically 
to  give 


V°>  "  Bv  {1  -  i/d0(y)  +  I1(y))}  (12) 

V3  TT 

where  y  =  »  3  K^R  and  Iq  and  are  modified  Bessel  Functions. 
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In  this  way  the  integral  in  equation  (II)  can  be  evaluated,  and 
the  equation  used  to  determine  E  for  any  given  cylinder  o*  radius 
R  and  temperature  T  The  current  I  is  then  obtained  from  the 
conductivity  equation  (13) 


I/F,  >■  ttR^o 

Curves  computed  using  equations  (11),  (12)  and  (13)  are  given  in 
Figs.  5,  6  and  7  for  pressures  of  2.5,  4  and  11  atmospheres 
respectively.  The  agreement  of  the  central  temperature  and  radius 
predicted  for  pressures  of  4  and  U  atmospheres  for  E  =  29  and 
I  -  1140  with  the  temperature  profiles  of  Figs.  1  and  2  is  within  5% 
Fcr  p  -  2.5  Atmospheres,  the  point  at  E  =  29  volt  cm”1 
and  I  -  H40  amperes  appears  to  be  .lust  in  the  folded  region  of 
the  curves  of  Fig.  5  where  multiple  solutions  are  possible.  The 
existence  of  multiple  solutions  has  been  investigated  cut  for 
other  input  material  functions  and  for  a  different  gas  and  will  be 
reported  later.  For  xenon  at  2.5  Atmospheres,  for  E  greater  than 

^0  and  1  *  1140  the  only  solution  that  exists  is  for  a  conduction 
dominated  arc  where  the  central  temperature  is  controlled  by 
thermal  conduction  and  not  by  radiation.  At  E  =  29  volt  cm"1 
I  -  1140  amperes  the  only  solution  that  was  found  was  that  of  a 
conduction  dominated  arc  w^th  a  central  temperature  above  19,000°K. 

At  these  temperatures  the  input  material  functions  are  unknown. 

No  detailed  search  was  made  for  the  existence  of  a  solution  corres¬ 
ponding  to  a  radiation  dominated  arc  for  this  value  of  E  and  I. 
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III.  COMPARISON  WITH  EXPERIMENT 
1.  Steady  state  condition 

For  a  tube  of  R  ■  0.63  cm  the  peak  values  of  I,  E  and 

temperature  that  ha^e  been  measured  are  E  ■  29  volt  cm  \  I  *  1140 

amps,  central  temperature  *  11,500°K  and  pressure  «*  4.2  atmosphere .  ^ 

3T 

As  these  quantities  apply  for  a  time,  corresponding  to  —  ■  0 

O  t 

it  »hould  be  appropriate  to  compare  with  the  steady  state  theoretical 
solution.  Upper  and  lower  limits  of  the  possible  gas  pressure  are 
11  and  2.5  Atmospheres. 

Results  obtained  from  Fig.  6,  and  for  similf  r  curves  at 
4  Atmospheres  but  with  different  values  of  a  and  are  given 
in  Table  2. 


Table  2  (p  -  4,  E  -  29,  I  -  1140) 

Input  Output 


0 

K 

V 

R 

T 

De  Votto 

Values 

Normal  Continuum 

0.4 

13,300 

De  Votto 

Values 

2  x  Continuum 

0.42 

11,700 

1/2  x  De 

Votto 

Normal  ConMnuum 

0.6 

11,700 

1/2  x  De 

Votto 

2  x  Continuum 

0.62 

10,800 
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If  lines  are  included  in  the  spectral  absorbtivity  data, 
the  integrated  emission  coefficient  may  be  increased  by  a  factor  of 
two.  i  jwever,  as  Is  seen  from  Table  2,  the  computed  radius  is  stij.1 
only  0.42  cm  and  a  would  stillneed  to  be  decreased  for  the  arc 
r.o  fill  the  tube. 

Unfortunately  measurements  of  the  discharge  radius  are 
difficult  to  make.  Furthemore.  no  accurate  experimental  measurements 
of  o  for  xenon  are  known.  !«.  is  seen  that  the  theoretical  predic¬ 
tions  are  in  very  approximate  agreement  vith  experiment,  but  if  the 
discharge  does  indeed  fill  the  discharge  tube  it  would  appear  that 
the  De  Vott  >  values  of  u  are  too  high.  However  an  analysis  of  the 
time  dtoendence  of  the  experimental  parameters  suggests  tha  the 
discharge  does  not  fill  the  tube.  Therefore  it  should  not  be 
concluded  that  ue  believe  that  the  a  or  :  data  is  in  error. 

2.  Time  dependence  of  the  discharge 

4 

From  the  time  dependence  of  the  experiment ally  observed 

O 

values  of  E,  I  and  radiation  intensity  at  8232A,  the  curves  of 
Fig.  ?  and  9  have  been  prepared.  In  Fig.  8  it  is  s°enthat  the 
turves  for  differing  discharge  energies  are  in  reasonable  coincidence 
for  the  latter  portion  of  the  discharge  time.  However  it  is  seen 
in  Ftg.  9  that  tne  radiation  intensity  at  a  particular  value  of  I 
varies  markedly  for  disiha^ges  of  differing  input  power.  If  the 
discharge  were  to  be  wall  stabilized  and  effectively  in  the  steady 
state  condition,  the  value  of  E  and  central  temperature  (and  ’•.hus 
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the  radiation  intensity)  would  be  determined  for  any  >ralue  of 
R  and  I.  If  the  experimeital  measurements  are  correct,  it  appears 
that  the  arc  cannot  be  regarded  as  being  wall  stabilized  and  in 
the  steady  state  at  times  after  t  e  peak  radiation  intensity. 

It  is  known  that  the  pressure  inside  of  the  discharge 
tube  is  determined  not  only  by  the  core  temperature,  but  is  con¬ 
trolled  to  a  large  extent  by  the  temperature  of  the  gas  near  the 
tube  wall.  It  is  suggested  that  this  gas  is  continually  heated 
throughout  the  duration  of  the  pulse,  thus  increasing  the  pressure. 
The  centra.1  temperature,  hich  is  determined  from  equation  (11), 
is  thus  lower  at  a  given  I  for  pulses  of  high  power  because  the 
gas  pressure  is  higher.  In  this  way  it  may  be  possible  to  explcln 
the  radiation  measurements  of  Fig.  9.  It  would  then  be  expected 
that  the  radius  of  the  central  portion  of  the  arc  does  not  fill 
the  tube,  and  increases  throughout  the  duration  of  the  pulse. 

Measurements  of  the  temperature  profile  at  peak  power 
are  consistent  with  the  theoretical  predictions  that  the  tempera¬ 
ture  profile  is  very  flat. 


IV.  CONCLUSIONS 

For  the  experimental  parameters  typical  for  the  operation  of  xenon 
flash  tubes,  the  following  generalizations  can  be  made. 

1.  The  temperature  trofiles  are  very  flat  and  are  radiation 
dominated,  provided  the  gas  pressure  equals  or  is  greater 
than  4  Atmospherea. 
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2.  Except  for  very  high  powers  where  the  central  temperature  is 
above  14,000°K,  the  effect  cf  radiation  of  wavelengths  less 
than  100UA  can  be  neglected. 

3.  Self  absorption  effects  in  the  visible  continuum  generally 
effect  the  central  temperature  by  less  than  1000°K. 

4.  The  constant  temperature  model  in  which  account  is  taken  of 
absorption,  can  be  used  to  predict  the  central  temperature 
and  outer  radius  of  wall  stabilized  arcs  to  within  10%. 

5.  Of  the  three  methods  used  to  compute  temperature  profiles  of 
steady  state,  wall  stabilized  arcs  from  input  material  functions, 
it  appears  that  the  time  dependent  method  using  an  exact 
radiation  intensity  calculation  is  the  most  economic, 

6.  Experimental  measurements  indicate  that  the  arc  discharge  is  not 
wall  stabilized  and  that  time  dependent  studies  are  needed. 

V.  POSSIBLE  FUTURE  WORK 
1.  Experimental 

(1)  Time  dependent  studies  of  the  radiation  intensity  as  a 
function  of  radius  are  needed  to  determine  whether  the 
radius  varies  with  time. 

(2)  Improved  measurements  of  the  pressure  as  a  function  of  time. 

(3)  Measurements  of  temperature  profiles  of  a  cascade  arc  in 
xenon.  Such  a  measurement  would  enable  a  determination 

of  electrical  conductivity,  which  is  a  crucial  input  parameter 
in  the  theoretical  calculation  of  temperature  profiles- 
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2.  Theoretical 

(1)  If  the  pressure  and  E  were  known  as  a  function  of  time,  a 
time  dependent  study  using  equation  (10)  could  be  made  to 
predict  the  current  and  the  temperature  profile  as  i 
function  of  time.  The  self-absorption  integral  of 
equation  (10)  could  be  neglected  giving  very  considerable 
imnrovemenf  in  computation  time  and  an  error  of  only  a 
few  hundred  degrees  in  the  central  temperature.  It  is 
proposed  that  for  a  given  I  the  temperature  of  the  gas 
near  the  quartz  walls  which  contain  the  arc,  is  dependent 
on  the  total  power  fed  into  the  pulse,  but  tKat  the 
temperature  profile  of  the  arc  core  at  late  times  is 
relatively  insensitive  to  the  initial  temperature  profile 
that  is  assumed  for  the  calculations.  The  gas  pressure 
as  a  function  of  time  could  be  derived  from  the  computed 
temperature  profiles  and  compared  with  the  input  pressure 
as  a  function  of  time. 

(2)  A  more  involved  computation  could  be  carried  out  solving 
equations  (14),  (15)  and  (16) 

C  p  (ff  +  v.  V  T)  -  V.(kVT)  +  v  rr  +  v  v.Vv 
p  ot  dt  — 

-  oE^  -  4tt  -  —  =  0  (14) 

p  dt 

3  v 

0  — +  P  V.v  v  -  -V  P  (15) 


||  +  V  .  pv  =  0 


(16) 
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As  p  and  p  determine  T;  p,  p,  v  and  T  can  be  determined  as  a 
function  of  time  and  compared  with  experiments. 
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in  the  appendices,  to  calculate  the  transport  properties,  the  spectral 
absorptivities  for  the  lines  and  the  continuum  of  xenon,  and  the  spectral 
radiance  and  temperature  profiles  in  cylindrical  arcs.  Also  included  as 
an  appendix  is  a  theoretical  analysis  of  the  xenon  arc  using  radiative 
transport  techniques  developed  in  other  studies. 
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